Tumor Radiomic Features on Pretreatment MRI to Predict Response to Lenvatinib plus an Anti-PD-1 Antibody in Advanced Hepatocellular Carcinoma: A Multicenter Study

Hepatocellular carcinoma (HCC) is one of the most prevalent malignant tumors and one of the leading causes of cancer mortality in the world and China [1, 2]. Most patients with HCC are diagnosed at an advanced stage, non-amenable to curative treatment [3].

In recent years, combination strategies including targeted therapy plus immunotherapy, such as atezolizumab plus bevacizumab and sintilimab plus a bevacizumab biosimilar [4], have been approved for the first-line treatment of advanced HCC (the latter only in China). Furthermore, combination therapy with lenvatinib plus pembrolizumab for the first-line treatment in advanced HCC was investigated in a phase Ib trial (KEYNOTE-524 study) and showed a promising objective response rate (ORR) of 36% (according to Response Evaluation Criteria in Solid Tumors, version 1.1 (RECIST v1.1)), with promising progression-free survival (PFS) and overall survival (OS) [5]. Similar ORRs have been reported for lenvatinib in combination with a range of different anti-programmed cell death protein 1 (PD-1) antibodies [6]. Although LEAP-002 study (lenvatinib plus pembrolizumab vs. lenvatinib as first-line therapy for advanced HCC, NCT03713593) is a negative trial in terms of PFS and OS, it is also noted that ORR was much higher in the combination arm than the lenvatinib monotherapy arm [7], which implies there may be a role of this combination therapy in neoadjuvant therapy. Indeed, preliminary investigation suggested a great value of this combination treatment in either conversion therapy or neoadjuvant therapy settings [8‒10]. A multicenter prospective clinical trial has been initiated in China to further investigate the efficacy of this combination in the neoadjuvant setting (NCT05389527). Therefore, it is also important to predict efficacy of this combination therapy, especially in the neoadjuvant setting when we already have multiple choices for systemic therapies.

Although the combination of lenvatinib plus an anti-PD-1 antibody provides a relatively high ORR, the clinical benefit is limited to a subset of patients. As the incidence of grade 3/4 adverse events is ∼50% with this combination [5, 11], identifying patients who are more likely to respond before initiating treatment is of clinical significance. However, currently, no validated biomarker has been identified to predict the response to combination therapy with targeted therapy plus immunotherapy [12].

Biomarkers such as programmed death-ligand 1 (PD-L1) expression, tumor mutational burden, microsatellite instability [13], and mismatch repair deficiency [14] have been proposed to predict response to immunotherapy in HCC and other cancers. However, whether these biomarkers could predict the efficacy of combination treatment with targeted therapy plus immunotherapy remains unclear [12]. Furthermore, these are tissue-based biomarkers that require invasive biopsies or procedures to obtain tumor tissue samples and do not adequately characterize the spatial and temporal heterogeneity of tumors.

Radiomics is an emerging field that converts standard-of-care medical imaging into high-throughput, mineable, and quantitative features using a variety of predefined image-characterization algorithms [15, 16]. Medical images contain information that directly reflects the overall tumor burden and each tumor lesion [16]. Therefore, radiomics, which can be perceived as a “digital biopsy” of the entire tumor, provides a noninvasive, repeatable, and comprehensive view of tumor biology and heterogeneity, without the need for additional blood or tissue samples.

Radiomics have been shown to improve diagnostic [17], prognostic [18, 19], and predictive accuracy [20, 21] in HCC. Radiomics can also be used to predict response to transcatheter arterial chemoembolization (TACE) in HCC [22], neoadjuvant therapy in rectal and breast cancer [23, 24], systemic therapy in lung cancer [25], and to immune checkpoint inhibitors in solid tumors [26, 27]. However, to the best of our knowledge, no study has investigated the ability of radiomics to predict response to lenvatinib plus an anti-PD-1 antibody in advanced HCC.

Magnetic resonance imaging (MRI) provides a higher soft-tissue contrast than computed tomography (CT) images, and contrast-enhanced MRI remains the most reliable radiologic technique to evaluate response to systemic therapy in HCC [28]. Therefore, we hypothesized that a radiomic model to predict response to lenvatinib plus an anti-PD-1 antibody in unresectable or advanced HCC based on pretreatment contrast-enhanced MRI would have similar or incremental value to a clinicopathologic model. In this study, we aimed to develop a radiomic model based on pretreatment contrast-enhanced MRI to predict objective response to the first-line combination therapy with lenvatinib plus an anti-PD-1 antibody in patients with unresectable or advanced HCC, and validate this radiomic model using a multicenter dataset.

This was a retrospective, multicenter study involving 9 tertiary referral centers in China. Consecutive patients diagnosed with HCC according to local and international guidelines [29, 30] were eligible for inclusion if they received combination therapy with lenvatinib plus an anti-PD-1 antibody (see online suppl. methods for details of treatment; see online suppl. Result S1, online suppl. Table S9 and online suppl. Table S10 for details of treatment-related adverse events; for all online suppl. material, see www.karger.com/doi/10.1159/000528034) as a first-line systemic treatment between October 2018 and February 2022, had undergone a pretreatment contrast-enhanced MRI within 2 weeks before initiating combination therapy, and had tumor response assessments every 2 months (±2 weeks) via CT or MRI according to RECIST v1.1. Patients were also required to have an interval of at least 2 months between any previous therapy (e.g., TACE, radiofrequency ablation, hepatectomy) and the initiation of combination therapy, and have at least one tumor response assessment after initiating combination therapy. Different anti-PD-1 antibodies were permitted across the study population to reflect real-world practice more accurately and showed similar ORRs when combined with lenvatinib [6]. Exclusion criteria were incomplete clinicopathologic data; a history of cancer other than HCC; a gadoxetic acid-enhanced pretreatment MRI; inadequate MRI quality; intrahepatic tumor lesions that could not be measured or segmented in magnetic resonance (MR) images; and additional antitumor treatment after the initiation of combination therapy.

Eligible patients formed two cohorts: a training cohort recruited from an ongoing observational, prospective cohort study (NCT04639284) at one of the centers (Zhongshan Hospital, Fudan University), whose data were used to build the model; and an independent validation cohort recruited from the remaining 8 centers (online suppl. methods).

The objective of this study was to predict objective response (defined as a complete response or partial response as the best overall response) to the first-line combination therapy with lenvatinib plus an anti-PD-1 antibody with an area under the receiver operating characteristic (ROC) curve (AUC) of ≥0.8. To avoid model overfitting, the rule of thumb is that the number of predictors should remain within 1/20-1/6 of the sample size in the training cohort used to build a model [31]. A sample size of at least 33 patients (10 patients with objective responses and 23 patients without objective responses) in each cohort was required based on the following assumption: a power of 0.8, two-sided α of 0.05, alternative hypothesis of an AUC of 0.8 compared with the null hypothesis of an AUC of 0.5, and an expected ORR of ∼30% based on a previous report [6]. Sample size was calculated using PASS 2021 (NCSS, LLC, Kaysville, UT, USA).

MRI data acquisition and scan parameters are described in online supplementary methods and online supplementary Table S1, respectively.

Tumor segmentation was performed using 3D slicer (version 4.11). The regions of interest for all intrahepatic target tumors were manually drawn along the tumor boundary on arterial phase and delayed phase images by a radiologist (S-Y Dong), and 60 randomly selected tumors were re-segmented by a senior radiologist (S-X Rao) to test the robustness and reproducibility of extracted radiomic features. Both radiologists were blinded to clinical, laboratory, and response assessment results for all patients. A representative case of tumor contouring on MR images is shown in online supplementary Figure S1.

Image preprocessing and radiomic feature extraction were performed using PyRadiomics (v3.0.1) [32]. Briefly, MRI signal intensity was normalized followed by resampling and interpolation of voxels to minimize the effects of different MRI acquisition parameters and scanners, which are detailed in the online supplementary methods. Finally, 1,118 features were extracted for each of the arterial and delayed phases (2,236 in total). The extracted features in each image type are listed in online supplementary Table S2.

Radiomic features with an interobserver intraclass correlation coefficient (ICC) <0.8 were excluded to ensure robustness and reproducibility as previously reported [31] and were normalized using z-score method.

The framework for the development and validation of the predictive models is shown in Figure 1a.

Clinicopathologic features with a p value <0.2 in univariate logistic regression analyses along with features that may have clinical significance [33] were selected to build a clinicopathologic model. Radiomic features with an interobserver ICC ≥0.8 and a p value <0.1 in Student’s t test were selected by least absolute shrinkage and selection operator logistic regression with stratified 5-fold cross-validation, maximizing the AUC, to build a radiomic model. The features in the clinicopathologic model and the radiomic score (rad-score) calculated using the radiomic model was used to build the clinicopathologic-radiomic model.

All prediction models were built using a neural network model in the training cohort. To avoid overfitting and increase generalizability, hyperparameter optimization was performed. The neural network model and hyperparameter optimization are described in the online supplementary methods. The prediction models were validated in the independent validation cohort to assess generalizability following development in the training cohort.

Bootstrap resampling (n = 1,000) was used to calculate an AUC with 95% confidence intervals (CI) and compare the difference between two ROCs. The optimal cut-off threshold of prediction models was determined using ROC by maximizing the Youden index in the training cohort. Net reclassification index (NRI) was used to quantify how well a new model reclassifies patients compared with an old model, at their respective cut-off values, which is particularly useful if there is no significant difference between the AUCs for the new and old models [34]. A NRI >0 indicates the predictive ability of the new model is better than that of old model. Calibration curves were plotted to assess the calibration of models. Decision curve analysis was conducted to determine the clinical utility of the prediction models by quantifying the net benefits at different threshold values [35].

Patients were followed every 60 days (±7 days) after initiation of combination therapy. OS was calculated from the date of first dose of drug to death from any cause or censored on the last follow-up. PFS was calculated from the date of first dose of drug to the first documented disease progression or death from any cause.

Categorical variables were expressed as counts and percentages, and were compared using Pearson’s χ² analysis, Fisher’s exact test, or Mann-Whitney U test, as appropriate. Continuous variables were expressed as mean (± standard deviation) or median (interquartile range (IQR)) and were compared using Student’s t test or Mann-Whitney U test, as appropriate. The multivariate logistic regression was used to identify independent predictors. Survival curves were calculated using the Kaplan-Meier method and compared using the log-rank test. A p value <0.05 was considered statistically significant. All statistical analyses were performed using the R software (version 4.1.2; packages used listed in online suppl. Table S3).

In total, 170 eligible patients were enrolled, 124 in the training cohort and 46 in the independent validation cohort (Fig. 1b). Patients in the validation cohort had better Eastern Cooperative Oncology Group performance status than patients in the training cohort (p = 0.010); there were no significant differences in other baseline demographic and disease characteristics between the training and validation cohorts (Table 1). The ORRs were 29.8% (37/124) and 28.3% (13/46) in the training and validation cohorts, respectively (p = 0.991).

The baseline demographic and disease characteristics of patients with and without an objective response in the training and validation cohorts are summarized in online supplementary Tables S4 and S5. In the training cohort, a smaller proportion of patients who had an objective response had a hepatitis B virus DNA copy number >10³/mL compared with those who did not have an objective response (27.0% vs. 48.3%, p = 0.046). Patients who had an objective response also had lower baseline α-fetoprotein (AFP, ng/mL) levels than those who did not have an objective response (median (IQR), 80.4 (5.5–10,385) vs. 1,104 (74.85–14,559), p = 0.029). In the validation cohort, patients with an objective response had a higher proportion of an Eastern Cooperative Oncology Group performance status of 0 (100% vs. 66.7%, p = 0.020) and a lower proportion of AFP >400 ng/mL (23.1% vs. 60.6%, p = 0.049) than those without an objective response.

Clinicopathologic features with a p value <0.2 in univariate logistic regression analyses in the training cohort were hepatitis B virus DNA (>10³/mL vs. ≤10³/mL), extrahepatic disease (yes vs. no), and AFP (>400 ng/mL vs. ≤400 ng/mL) (Table 2). The additional features that may have clinical significance in predicting objective response were macrovascular invasion (yes vs. no), and sum of diameter of baseline intrahepatic target lesions in cm. The continuous form of AFP was used instead as the continuous form can provide more information than the categorical form.

The clinicopathologic model was built with 1, 4, and 2 neurons in hidden layers 1, 2, and 3, respectively. The diagram and formula of this model are shown in online supplementary Figure S2. The AUC (95% CI) using this model was 0.748 (0.656–0.840) and 0.702 (0.547–0.884) in the training and validation cohorts, respectively, after bootstrap resampling (n = 1,000) (online suppl. Fig. S3). The optimal cut-off threshold for this model was 0.614; accuracy, sensitivity, specificity, positive predictive value (PPV), and negative predictive value (NPV) for this model are summarized in Table 3. These results indicated that clinicopathologic features modestly predict the objective response to lenvatinib plus an anti-PD-1 antibody. Furthermore, the training and validation cohorts did not meet the study objective in the clinicopathologic model.

The heatmap expressions of extracted radiomic features in each patient are shown in Figure 2a, b. The details of radiomic feature selection and selected features are described in Figure 2c and online supplementary Result S2, respectively. Ultimately, 14 radiomic features selected in the arterial and 3 in the delayed phase were used to build the radiomic model.

The radiomic model was built with 2 neurons in the single hidden layer. The diagram and formula of this model are shown in online supplementary Figure S4. The AUC (95% CI) using this model was 0.886 (0.815–0.957) and 0.820 (0.648–0.984) in the training and validation cohorts, respectively, after bootstrap resampling (n = 1,000) (Fig. 3a). The rad-score was identified as an independent predictor for objective response in both the training and validation cohorts (online suppl. Table S6).

Rad-scores were significantly associated with objective tumor response both in the training and validation cohorts (p < 0.001 for both, Fig. 3b; the distributions of rad-scores in the training and validation cohorts are shown in Fig. 3c, d, respectively). The optimal cut-off threshold of rad-score was 0.504, and accuracy, sensitivity, specificity, PPV, and NPV for the radiomic model are summarized in Table 3. Patients with a rad-score >0.504 were predicted to have an objective response; otherwise, patients were predicted not to have an objective response. This cut-off discriminated between patients with and without an objective response in both the training and validation cohorts (p < 0.001 for both, Table 4).

There were no statistically significant differences in the clinicopathologic features between patients predicted to have an objective response and those predicted not to in the training cohort (online suppl. Table S7). Compared with patients predicted not to have an objective response in the validation cohort, a smaller proportion of patients predicted to have an objective response were AFP >400 ng/mL (11.1% vs. 59.5%, p = 0.022), and they had an earlier CNLC stage (p = 0.011) (online suppl. Table S8). The calibration curves showed good agreement between the outcomes predicted by the radiomic model and the actual clinical outcomes (online suppl. Fig. S5). The decision curves revealed that the radiomic prediction model can add net benefit than assuming all patients did or did not have an objective response in both the training and validation cohorts (Fig. 3e, f). These results demonstrated that radiomic features obtained from pretreatment MRI can predict objective response to lenvatinib plus an anti-PD-1 antibody with good discrimination, calibration, and clinical utility, independent on clinicopathologic features.

We incorporated rad-scores calculated using the radiomic model into the clinicopathologic model and built a clinicopathologic-radiomic model with 5, 3, and 2 neurons in the hidden layers 1, 2, and 3, respectively. The diagram and formula of this model are shown in online supplementary Figure S6. The AUC (95% CI) using this model was 0.987 (0.968–1.000) and 0.884 (0.762–1.000) in the training and validation cohorts, respectively, after bootstrap resampling (n = 1,000) (Fig. 3g).

The scores calculated using the clinicopathologic-radiomic model were significantly associated with objective tumor response in both the training and validation cohorts (p < 0.001 for both; Fig. 3h). The optimal cut-off for this model was 0.443, and accuracy, sensitivity, specificity, PPV, and NPV for this model are summarized in Table 3. Patients with a score >0.443 were predicted to have an objective response; otherwise, patients were predicted not to have an objective response. This cut-off threshold discriminated between patients with and without an objective response in both the training and validation cohorts (p < 0.001 for both, Table 4). The calibration curves showed good agreement between the outcomes predicted by clinicopathologic-radiomic model and the actual clinical outcomes (online suppl. Fig. S5). The decision curves revealed that this prediction model can add net benefit than assuming all patients did or did not have an objective response in both the training and validation cohorts (Fig. 3i, j).

Compared with the clinicopathologic model, the clinicopathologic-radiomic model had a higher AUC in the training cohort (0.987 vs. 0.748, p < 0.001) and a marginally higher AUC in the validation cohort (0.884 vs. 0.702, p = 0.074). According to NRI, the clinicopathologic-radiomic model improved prediction performance over the clinicopathologic model in both the training (NRI = 47.9%, p < 0.001) and validation cohorts (NRI = 41.5%, p = 0.025). These results indicated that radiomic features can provide incremental predictive value to clinicopathologic features for prediction of objective response to lenvatinib plus an anti-PD-1 antibody.

As of August 22, 2022, median follow-up was 15.2 (IQR: 8.5–22.2) months in the training cohort and 10.0 (IQR: 4.4–17.6) months in the validation cohort; median duration of treatment was 5.6 (IQR: 2.5–9.0) months in the training cohort and 6.0 (IQR: 2.6–14.8) months in the validation cohort. 81 (65.3%) of 124 patients in the training cohort and 28 (60.9%) of 46 patients in the validation cohort had progressive disease or died. The median OS was 20.1 (95% CI: 16.8–30.5) months in the training cohort and was not reached in the validation cohort; the median PFS was 10.7 (95% CI: 8.2–18.0) months and 11.4 (95% CI: 9.2–16.0) months in the training and validation cohorts, respectively (online suppl. Fig. S7).

Patients with a rad-score >0.504 were associated with a significantly longer median OS or median PFS than those with a rad-score ≤0.504 both in the training cohort (median OS: 30.5 (95% CI: 22.2–not evaluable (NE)) months versus 16.8 (95% CI: 14.3–22.1) months, hazard ratio (HR) (95% CI): 0.377 (0.192–0.740), p = 0.003; median PFS: 28.2 (95% CI: 12.6–NE) months versus 8.1 (95% CI: 6.1–10.9) months, HR (95% CI): 0.368 (0.203–0.667), p < 0.001) and validation cohort (median OS: NE (95% CI: NE–NE) months versus 14.4 (95% CI: 9.7–NE) months, HR (95% CI): 0.066 (0.001–0.504), p = 0.008; median PFS: NE (95% CI: 14.3–NE) months versus 10.5 (95% CI: 7.0–14.4) months, HR (95% CI): 0.177 (0.042–0.754), p = 0.009; Fig. 4).

Previous studies showed that modified albumin-bilirubin (mALBI) grade [36, 37] and neutrophil-to-lymphocyte ratio (NLR) [38] were associated with OS and PFS in lenvatinib treatment. In the drug combination setting of this study, mALBI grade was associated with OS in the training cohort (p = 0.018) but not in the validation cohort (p = 0.360); mALBI grade was not associated with PFS both in the training (p = 0.570) and validation (p = 0.270) cohorts (online suppl. Fig. S8).

The AUC (95% CI) using NLR to predict objective response was 0.514 (0.404–0.623) and 0.629 (0.400–0.826) in the training and validation cohorts, respectively, after bootstrap resampling (n = 1,000). The optimal cut-off threshold of NLR was 5.861 by Youden index using ROC. NLR was not associated with OS or PFS both in the training and validation cohorts (online suppl. Fig. S9). When using 2.548 (the median value of NLR) as the optimal cut-off threshold, NLR was not associated with OS in the training cohort (p = 0.095) but in the validation cohort (p = 0.014); NLR was not associated with PFS both in the training (p = 0.440) and validation (p = 0.350) cohorts (online suppl. Fig. S10). When using 3 [38] as the optimal cut-off threshold, NLR was associated with OS both in the training (HR (95% CI): 0.568 (0.347–0.929), p = 0.022) and validation (HR (95% CI): 0.291 (0.102–0.824), p = 0.014) cohorts, but not associated with PFS both in the training (p = 0.470) and validation (p = 0.380) cohorts (online suppl. Fig. S11).

In this study, we developed a radiomic model based on pretreatment MRI to predict objective responses to combination therapy with lenvatinib plus an anti-PD-1 antibody in patients with unresectable or advanced HCC, which was validated in a multicenter dataset, and indicated that radiomic features were associated with OS and PFS after initiating this combination therapy. To the best of our knowledge, this is the first study to report the radiomic analysis based on pretreatment contrast-enhanced MRI to predict response to combination therapy with lenvatinib plus an anti-PD-1 antibody in advanced HCC patients.

Tumor radiomic features were associated with response to anti-PD-1 antibodies or lenvatinib monotherapy plus TACE in advanced HCC in two recent studies [39, 40], providing important data to support the generation of the hypothesis of the present study. The first study proposed a radiomic model based on contrast-enhanced CT with a relatively small sample size (N = 58) from one center to predict response to anti-PD-1 antibodies monotherapy [39]. MRI has several advantages over CT; it can provide superior contrast resolution, does not rely on ionizing radiation [41], and remains the preferred modality to assess response to systemic therapy in HCC [28]. In the second study, tumor radiomic features extracted from pretreatment MRI predicted disease progression after lenvatinib monotherapy plus TACE [40]. It was also a single-center study with a small sample size (N = 61); the AUC was 0.71, which had a modest discrimination power.

All (17/17) the final radiomic features selected to predict objective response were from wavelet and Laplacian of Gaussian filtered images instead of original images, indicating that radiomics can be used to identify details and extract features on MRI that cannot be captured or quantified by the naked eye, and thus accurately reflect the biology of HCC. These high-dimensional features cannot be detected in original images or by the naked eye, but hold more detailed information about the tumors and more sensitively predicted treatment response than clinicopathologic features [24]. For example, the “3D_glszm_GrayLevelVariance” feature represents the discrete degree of gray level in tumor regions, which was associated with tumor heterogeneity (e.g., tumor cellularity, micro-necrosis, and inflammation). Previous studies have shown that tumor heterogeneity is related to tumor response, which is consistent with the results of this study [42, 43]. Furthermore, the “firstorder_Kurtosis” feature was associated the degree of tumor enhancement in the arterial phase, which reflects tumor vasculature. A previous study has shown that the enhancement degree of tumor in arterial phase before immunotherapy is related to its progression [44]. The remaining texture features may reflect the immune microenvironment of tumors, but need to be further studied with genomic or histological correlative data [31].

Consistent with the present study, the model based on clinicopathologic features alone did not satisfactorily predict objective response; however, the radiomic model robustly predicted objective response. Furthermore, incorporating radiomic features into the clinicopathologic model provided incremental predictive value for objective response, suggesting that clinicopathologic features and radiomic features may reflect different characteristics of tumors related to the response to combination therapy. The clinicopathologic-radiomic model had a marginally higher AUC than the clinicopathologic model (0.884 vs. 0.702, p = 0.074) in the validation cohort, which may be attributed to an insufficient sample size. However, using NRI, we confirmed the incremental predictive value of radiomic features compared with clinicopathologic features in predicting response. Tumor radiomic features may capture tumor biology and heterogeneity; however, the objective response is not only associated with tumor biology but also related to tumor burden and condition of patients. For example, according to the updated IMbrave150 data, patients with BCLC stage B disease had a higher objective response than all patients (43.0% vs. 29.8%) [45, 46]. Therefore, the clinicopathologic-radiomic model provided a modestly higher discriminative ability to predict the objective response than the radiomic model in each cohort.

As the prognosis of HCC patients depends on tumor burden and hepatic reserve function, previous studies demonstrated that the mALBI grade was a negative prognostic factor for OS in lenvatinib treatment against advanced HCC [36, 37]. However, mALBI grade was only associated with OS in the training cohort, but not in validation cohort, which might be attributed to limited patient number and immature survival data in the validation cohort, or different treatment regimen. Another study first indicated NLR <3 was a favorable factor for OS in HCC patients treated with lenvatinib [38]. In this study, NLR <3 was demonstrated to be associated with OS both in the training and validation cohorts for lenvatinib plus an anti-PD-1 antibody, but not associated with PFS.

This study has several limitations. First, it had a limited sample size with potential selection bias. The sample size requirement for developing the radiomic model in the training cohort was met, but more data are required to optimize and improve models. Second, the majority of the study population had hepatitis B virus-related HCC. As heterogeneity exists in patients with HCC of different etiologies [47], the proposed radiomic model may not be applied to HCC patients with nonalcoholic fatty liver disease or infected with hepatitis C virus. Additional studies are required to compare the differences in radiomic features among patients with HCC of different etiologies. Third, due to the retrospective nature of this study, the MR image acquisition parameters were not standardized across different hospitals, although this study demonstrated that MR image preprocessing could overcome this influence; furthermore, tumor tissues were not prospectively collected to explore the biological meaning of radiomic features. Fourth, manual segmentation of tumors is relatively subjective, especially for tumors with blurry edges, but ICC was used to reduce subjectivity and ensure robustness and reproducibility as previously reported [31]. Automatic or semiautomatic segmentation tools may reduce the interobserver ICC. Furthermore, only radiomic features from arterial and delayed phases were analyzed. Although characteristics of HCC on imaging are mainly reflected on these two phases, radiomic features from other sequences in MRI should be analyzed in the future. Fifth, different anti-PD-1 antibodies were used across the study population, which could increase internal variability with a limited sample size; however, lenvatinib in combination with different anti-PD-1 antibodies were reported to have similar ORRs [6] and could reflect real-world practice.

In conclusion, this study demonstrated that tumor radiomic features extracted from pretreatment MRI can be used to predict objective response to combination therapy with lenvatinib plus an anti-PD-1 antibody in patients with unresectable or advanced HCC, and provide incremental predictive value over clinicopathologic features, and are associated with OS and PFS after initiation of this combination regimen. A prospective study with a unified protocol to investigate the performance of radiomic and clinicopathologic-radiomic models and to explore their usage in patients receiving existing first-line systemic treatments for advanced HCC is warranted.

We thank the patients and their families. We also thank the China liver cancer study group young investigators (CLEAP) for their joint efforts for this work.

The study protocol was conducted in accordance with the principles of the World Medical Association Declaration of Helsinki and was approved by the Zhongshan Hospital Research Ethics Committee (Approval Numbers: B2021-210). The requirement for informed consent was waived (retrospective computational analysis of images) by the Zhongshan Hospital Research Ethics Committee.

Hui-Chuan Sun has received honorarium or lecture fees from Roche, Bayer, MSD, Eisai, Hengrui, Innovent, TopAlliance, Abbott, Beigene, Gilead, and Zelgen during the last 5 years. Jia Fan is an Editorial Board Member of Liver Cancer. All remaining authors have no conflicts of interest to declare.

This work was supported by the Leading Investigator Program of the Shanghai municipal government (17XD1401100), the National Key Basic Research Program (973 Program, 2015CB554005) from the Ministry of Science and Technology of China, the National Natural Science Foundation of China (81871928), the Special Research Fund for Liver Cancer Diagnosis and Treatment from the China Anti-Cancer Association (H2020-008), and the Clinical Research Special Fund of Zhongshan Hospital, Fudan University (2020ZSLC71) to Hui-Chuan Sun.

Study concept and design, drafting of the manuscript: Bin Xu and Hui-Chuan Sun; acquisition of data: Bin Xu, San-Yuan Dong, Xue-Li Bai, Tian-Qiang Song, Bo-Heng Zhang, Le-Du Zhou, Yong-Jun Chen, Zhi-Ming Zeng, Kui Wang, Hai-Tao Zhao, Na Lu, Wei Zhang, Xu-Bin Li, Su-Su Zheng, Guo Long, Yu-Chen Yang, Hua-Sheng Huang, Lan-Qing Huang, Yun-Chao Wang, Xiao-Dong Zhu, Cheng Huang and Ying-Hao Shen; analysis and interpretation of data and critical revision of the manuscript for important intellectual content: all authors; statistical analysis: Bin Xu, San-Yuan Dong, Fei Liang and Hui-Chuan Sun; obtained funding: Hui-Chuan Sun; administrative, technical, or material support, study supervision: Xue-Li Bai, Tian-Qiang Song, Bo-Heng Zhang, Le-Du Zhou, Yong-Jun Chen, Zhi-Ming Zeng, Kui Wang, Hai-Tao Zhao, Jian Zhou, Meng-Su Zeng, Jia Fan, Sheng-Xiang Rao, and Hui-Chuan Sun.

Raw data of MRI used in this study are unable to be shared as the requirement of Hospital Research Ethics Committees. Derived data (e.g., extracted radiomic features) are available upon reasonable request from the corresponding authors.

Bin Xu and San-Yuan Dong contributed equally to this work.