Use of Radiomics to Predict Response to Immunotherapy of Malignant Tumors of the Digestive System.

BACKGROUND Despite the promising results of immunotherapy in cancer treatment, new response patterns, including pseudoprogression and hyperprogression, have been observed. Radiomics is the automated extraction of high-fidelity, high-dimensional imaging features from standard medical images, allowing comprehensive visualization and characterization of the tissue of interest and corresponding microenvironment. This study assessed whether radiomics can predict response to immunotherapy in patients with malignant tumors of the digestive system. MATERIAL AND METHODS Computed tomography (CT) images of patients with malignant tumors of the digestive system obtained at baseline and after immunotherapy were subjected to radiomics analyses. Radiomics features were extracted from each image. The formula of the screened features and the final predictive model were obtained using the Least Absolute Shrinkage and Selection Operator (LASSO) algorithm. RESULTS Imaging analysis was feasible in 87 patients, including 3 with pseudoprogression and 7 with hyperprogression. One hundred ten radiomics features were obtained before and after treatment, including 109 features of the target lesions and 1 of the aorta. Four models were constructed, with the model constructed from baseline and post-treatment CT features having the best classification performance, with a sensitivity, specificity, and AUC of 83.3%, 88.9%, and 0.806, respectively. CONCLUSIONS Radiomics can predict the response of patients with malignant tumors of the digestive system to immunotherapy and can supplement conventional evaluations of response.

Background

Immunotherapy has become a treatment of choice for patients with refractory or recurrent tumors. Unlike cytotoxic drugs, immunotherapeutic agents activate the anti-tumor activities of the immune system [1-4]. Immune checkpoint antibodies identified to date include antibodies to cytotoxic T lymphocyte-associated antigen-4 (CTLA-4) and programmed cell death-1/programmed cell death-ligand 1 (PD-1/PD-L1). These antibodies have shown efficacy in the treatment of patients with melanoma, non-small cell lung cancer, renal cancer, and head and neck squamous cell carcinoma, as well as in microsatellite instability high (MSI-H) and mismatch repair deficient (dMMR) solid tumors.

Immunotherapeutic agents, however, can affect the immune system and alter anti-tumor responses, giving rise to unique response patterns, including pseudoprogression [5] and hyperprogression [6,7]. Pseudoprogression is defined as a reduction in tumor burden after its increase or the appearance of new lesions [8,9], whereas hyperprogression is defined as a tumor growth rate >2-fold higher after than before immunotherapy [6,7]. These response patterns make the evaluation of immunotherapy more complicated. For sample, some patients with increased tumor burdens or new lesions may benefit from continued immunotherapy, whereas others experience progressive disease (PD), even hyperprogression, suggesting that they should be switched to alternative treatments. It is difficult to distinguish between patients who do and do not benefit from immunotherapy, including those patients who experience exacerbation of their condition while on immunotherapy, based on current response criteria.

The rapid development of medical informatics technology has given rise to radiomics, a method involving the automated extraction and analysis of a large number of high fidelity, advanced quantitative imaging features from medical images [10]. These medical images are used to comprehensively and non-invasively monitor tumor development and responses to treatment, providing a reliable solution to the problem of tumor heterogeneity [10,11]. Radiomics has been applied to the diagnosis and clinical staging of various types of tumors, including lung, colorectal, breast cancers, and gliomas, as well as to evaluating responses to treatment and prediction of patient prognosis [12]. Radiomics is especially important in patients receiving targeted therapy and immunotherapy, as necrosis of tumor cells may not occur during the course of treatment. Thus, evaluating response based on tumor diameter likely will not accurately reflect the therapeutic efficacy of these treatments. The present study assessed whether radiomics could predict response to immunotherapy in patients with malignant tumors of the digestive system.

Material and Methods

Data were retrospectively collected from patients with malignant tumors of the digestive system who had been treated with anti-PD-1/PD-L1 antibody, alone or in combination with anti-CTLA-4 antibody, in the Department of Gastrointestinal Oncology of Peking University Cancer Hospital from July 2016 to November 2017. All patients were evaluated by enhanced computed tomography (CT) scans of the chest, abdomen, and pelvis. The responses of tumors were evaluated by irRECIST. Newly diagnosed lesions (≤2 lesions per organ, ≤5 lesions in total) were included in the total measured tumor burden (TMTB).

Pseudoprogression was defined as irPD during evaluation, followed subsequently by irCR, irPR, or irSD. The dynamic growth of tumors was evaluated by tumor growth kinetics (TGK). TPRE, T0, and TPOST represent scanning times before baseline, at baseline, and at the first evaluation after baseline, respectively. TMTBPRE, TMTB0, and TMTBPOST represent TMTB before baseline, at baseline, and at the first evaluation after baseline, respectively. TGKPRE was defined as the change in TMTB per unit time before immunotherapy and was calculated as (TMTB0–TMTBPRE)/(T0–TPRE). Similarly, TGKPOST was defined as the change in TMTB per unit time after immunotherapy and was calculated as (TMTBPOST–TMTB0)/(TPOST–T0). The ratio of TGKPOST to TGKPRE was defined as the TGK ratio (TGKR). Hyperprogression was defined as a TGKR ≥2 or a ≥50% increase in TMPBPOST from baseline.

Imaging data were retrospectively collected through the picture archiving and communication system (PACS) of Peking University Cancer Hospital. All patients had undergone enhanced and unenhanced CT scanning of the chest and abdomen at baseline and after treatment. CT scanning was performed using a GE Discovery 750HD, GE Lightspeed VCT, or Philips iCT scanner, with a tube voltage of 120 kV, an automatic current, a tube speed of 0.8–1.0 s/r, and a collimation of 64*0.625 mm, with reconstruction layer thickness of 5 mm and layer spacing of 5 mm. Patients were injected with 80–100 ml of the contrast agent iohexol (300 mg/ml) through the cubital vein at a flow rate of 3.0 ml/s. Abdominal and pelvic enhanced CT scans were taken in the hepatic arterial phase at 25–30 s, the portal venous phase at 60–70 s, and the equilibrium phase at 90–100 s, whereas chest enhanced scans were taken 25–30 s after administration of contrast agent. Patients missing imaging data at the corresponding time points and those evaluated by other imaging methods (e.g., MRI) were excluded.

Regions of interest (ROIs) of target lesions in portal venous phase CT images were delineated along the outline of the target lesion at the maximum level using ITK-SNAP software. For reference, the aorta was delineated at the same time. Each ROI was manually segmented by an experienced radiographer. The digital images and sketched ROIs were subsequently uploaded to MATLAB software for feature extraction. The image features extracted in this study were divided into 3 categories with a total of 109 features, including 15 first-order gray-level statistical features, 16 geometric features, and 78 second-order texture features (Table 1).

Table 1

Classification of the extracted imaging features.

First-order gray-level statistical features (n=15)	Geometric features (n=16)	Second-order texture features (n=78)
Maximum gray valueMinimum gray valueMedianSumAverage valueStandard deviationVarianceKurtosisSkewness etc.	VolumeLong axis lengthShort axis lengthSurface areaEccentricityExtensionVolume of the cuboidMaximum perimeterDirection etc.	Gray-level co-occurrence matrix (glcm) (n=16)Gray-level run-length matrix (glrlm) (n=26)Gray-level size zone matrix (glszm) (n=26)Neighborhood gray-tone difference matrix (ngtdm) (n=10)

Feature screening involved 3 steps. In the first step, significant screening, t values between positive and negative samples were calculated for each feature using the equation:

where χ̆ indicates the mean, S2 indicates the sample variance, and n indicates the sample size. The t values were arranged from the largest to the smallest. The feature with the larger t value was selected, and the feature with the smaller t value was discarded.

In the second step, colinearity screening, the covariance between features was calculated by correlation coefficients, with colinearity defined as features with a correlation coefficient >0.9. Correlation coefficients were calculated using the formula:

where S indicates the standard deviation. If 2 features showed collinearity, the feature with the larger t value was retained and the feature with the smaller t value was rejected.

In the third step, collaborative screening features were screened and a final linear prediction model constructed using the Least Absolute Shrinkage and Selection Operator (LASSO) [13]. The LASSO algorithm suppresses the coefficients of some features by adjusting their weight parameters λ using the formula:

where βTχ̆ is the abbreviated form of wTχ+b, and the linear regression model predicts z=wTχ+b.

For model establishment and verification, the LASSO algorithm provided the linear combination formula of the screened features and the final predictive model. The predictive model was trained by a logistic regression model, followed by 5-fold cross-validation of the data set to verify the predictive ability of the model. Receiver operating characteristic (ROC) curves were constructed, and the performance of the model was assessed by determining its sensitivity, specificity, and the area under the ROC curve (AUC).

Results

During the study period, 112 patients with malignant tumors of the digestive system were treated with anti-PD-1/PD-L1 antibody; the characteristics of these patients are summarized in Table 2. First response was evaluated in 97 patients by irRECIST, with 24 achieving irPR, 26 having irSD, and 47 having irPD. No patient achieved irCR. Four patients were identified with pseudoprogression during follow-up, and 8 had hyperprogression, defined as a TGKR ≥2 and a >50% increase in TMTBPOST. The patient flowchart is shown in Figure 1.

Table 2

Demographic and clinical characteristics of patients treated with immunotherapy.

	Patients (N=112)
Sex
Male	69 (71.1%)
Female	28 (28.9%)
Location
Stomach	34 (35.1%)
Esophagus	21 (21.6%)
Colorectum	20 (20.6%)
Pancreas	12 (12.4%)
Hepatobiliary	8 (8.2%)
Intestine	2 (2.1%)
Histological type
Adenocarcinoma	51 (52.6%)
Squamous carcinoma	20 (20.6%)
Neuroendocrine carcinoma	20 (20.6%)
Hepatocellular carcinoma	4 (4.1%)
Cholangiocarcinoma	2 (2.1%)
Radical operation
No	44 (45.4%)
Yes	53 (54.6%)
Radiotherapy
No	69 (71.1%)
Yes	28 (28.9%)
Targeted therapy
No	82 (84.5%)
Yes	15 (15.5%)
ECOG performance status
0	48 (49.5%)
1	49 (50.5%)
PD-L1 status
Negative	21 (42.0%)
Positive	29 (58.0%)
MMR status
pMMR	37 (63.8%)
dMMR	21 (36.2%)
Immunotherapy type
PD-1	65 (67.0%)
PD-1+CTLA-4	3 (3.1%)
PD-L1	24 (24.7%)
PD-L1+CTLA-4	5 (5.2%)

ECOG – Eastern Cooperative Oncology Group; MMR – mismatch repair protein; pMMR – proficient MMR; dMMR – deficient MMR.

Figure 1

Patient flowchart of the selection process.

Of the 97 patients evaluated for first response, 10 were excluded because they were analyzed by other imaging methods (such as MRI) or there were no imaging data in the PACS. Of the 87 patients with feasible imaging analysis, 3 had pseudoprogression and 7 had hyperprogression. First responses included iPR in 17 patients, irSD in 24, and irPD, in 36. The patients were divided into 2 groups, those who benefitted from treatment (defined as those with pseudoprogression, irPR, and irSD) and those who did not (defined as those with hyperprogression and irPD).

The features of CT images at baseline and at first response evaluation were extracted and normalized relative to the mean signal of the aorta. A total of 110 features were obtained before and after treatment, including 109 features of the target lesions and 1 feature of the aortas.

Model 1 was established using the baseline CT features, as screened by the LASSO method (Figure 2A). Of the 110 features identified, 30 were selected (Table 3). The sensitivity, specificity, and AUC of model 1 were 64.3%, 39.5%, and 0.646, respectively (Table 4). Similarly, model 2 was established using the CT features screened at first evaluation. One of the 110 features was selected (Figure 2B). The sensitivity, specificity, and AUC of model 2 were 43.2%, 81.0%, and 0.750, respectively. Model 3 was established using the CT features at baseline and at first evaluation of response. Three of the 220 features were selected (Figure 2C). The sensitivity, specificity, and AUC of model 3 were 83.3%, 88.9%, and 0.806, respectively. Finally model 4 was established using the difference between CT features at baseline and at first evaluation of response. Ten of the 110 features were selected (Figure 2D). The sensitivity, specificity, and AUC of model 4 were 83.4%, 66.7%, and 0.806, respectively.

Figure 2

Screening of features using the LASSO method: (A) CT features at baseline; (B) CT features at first evaluation; (C) CT features at baseline and first evaluation; (D) difference between CT features at baseline and at first evaluation. The x-axis shows the parameter λ and the corresponding number of features (degrees of freedom), and the y-axis shows the average AUC values. The dotted line corresponds to λ at the maximum AUC and the number of features.

Table 3

Constructed models and their extracted CT features.

Model	First-order gray-level statistical features	Geometric features	Second-order texture features
1	Minimum gray valueMeanVarianceGray histogram features (n=3)	2D_eccentricity3D_surface area	GLCM (n=4)GLRLM (n=6)GLSZM (n=9)NGTDM (n=3)
2	–	–	GLCM (n=1)
3	–	3D_short axis length	GLCM (n=1)GLRLM (n=1)
4	–	3D_long axis length	–

Table 4

Statistical results of the constructed models.

Model	Number of features	AUC	Sensitivity (%)	Specificity (%)
1	30	0.646	64.3	39.5
2	1	0.750	43.2	81.0
3	3	0.806	83.3	88.9
4	1	0.806	83.4	66.7

Of these 4 models, model 1 had the lowest specificity and model 2 had the lowest sensitivity. Although the AUCs of models 3 and 4 were identical, 0.806, the specificity of model 4 was slightly lower. Therefore, of these 4 models, model 3 showed the best classification performance for predicting patient response to immunotherapy.

Only 3 patients with feasible imaging analysis had pseudoprogression, a number too small to construct a model. The findings in these 3 patients were therefore used to test the 4 models. Model 3 had an AUC of 0.736 and was correct in predicting responses in 2 of 3 patients. In comparison, model 4 had an AUC of 0.760 and was also correct in predicting responses in 2 of 3 patients. The patients who experienced tumor progression on treatment were further divided into 2 subgroups, one with hyperprogression and the other with irPD. A model was established using CT features at baseline and during initial evaluation of responses. Of the 220 features identified, 3 were selected (Table 5). The AUC for the prediction model was 0.877. A second model was formulated based on the difference between CT features at baseline and at first evaluation of response. One of the 110 features was selected (Table 5), and the AUC for the prediction model was 0.774. The features of the 2 models indicate that the maximum gray value has great significance for the determination of hyperprogression. This factor belongs to the first-order gray-level statistics feature, reflecting the distribution of individual pixel values without considering the spatial relationship.

Table 5

Models for hyperprogression and their extracted CT features.

Model	First-order gray-level statistical features	Geometric features	Second-order texture features
CT features at baseline and first evaluation	Maximum gray value(Baseline)Maximum gray value(Evaluation)	–	NGTDM (n=1)
Difference between CT features at baseline and first evaluation	Maximum gray value	–	–

The Kaplan-Meier survival curves generated using the 4 models are shown in Figure 3. Median overall survival was significantly lower in patients who progressed on treatment than in those who benefitted from treatment as shown using models 1 (8.7 vs. 25.5 months, p=0.005) (Figure 3A) and 2 (7.4 vs. 25.5 months (p=0.006) (Figure 3B). Models 3 (6.2 vs. 13.8 months; p<0.001) (Figure 3C) and 4 (6.2 vs. 13.8 months (p=0.001) (Figure 3D) also showed that median overall survival was significantly lower in patients who progressed than in those who benefitted from treatment.

Figure 3

Kaplan-Meier curves of patient overall survival based on the 4 models: (A) model 1, constructed from baseline CT features; (B) model 2, constructed from CT features at first evaluation; (C) model 3, constructed from CT features at baseline and at first evaluation; (D) model 4, constructed from the difference in CT features at baseline and at first evaluation.

Discussion

This study explored the usefulness of radiomics in evaluating the responses of tumors to immunotherapy. Four predictive models were constructed. Model 1 was constructed from CT features at baseline, whereas model 2 was constructed from CT features at first evaluation, with model 2 having fewer CT features and a higher AUC. The selected feature in the latter model was the heterogeneity of the gray-level of the co-occurrence matrix (2D_GLCM_Homogeneity), which has an unevenness of texture that reflects the efficacy of immunotherapy. However, its sensitivity was only 43.2%, indicating that the rate of missed diagnoses by model 2 was high. Model 3, constructed from CT features at baseline and at first evaluation, resulted in the extraction of 3 features. The sensitivity, specificity, and AUC of model 3 were 83.3%, 88.9%, and 0.806, respectively, with the sensitivity and AUC being superior to those of model 2. Model 4 was constructed from the difference between CT features at baseline and at first evaluation. One feature, variation in 3D long axis length, was selected. The AUC of this model was identical to that of model 3, but its specificity was slightly lower, 66.7%. Thus, of the 4 models tested, model 3 showed the optimal classification performance for predicting patient responses to immunotherapy.

To assess the prognostic value of these models, Kaplan-Meier survival curves were determined for each. Median overall survival differed significantly in patients who progressed on treatment and those who benefitted from treatment according to all 4 models. The difference between the 2 groups was largest in model 1 (8.7 vs. 25.5 months). Thus, the models constructed by radiomics could separate patients who did and did not benefit from immunotherapy.

Only 3 patients in this cohort had both pseudoprogression and feasible imaging analysis, a number too small to construct a model. Testing of these patients with models 3 and 4 yielded AUCs of 0.736 and 0.760, respectively. Of these 3 patients, 2 could be assigned to the group that benefitted from treatment. Testing of patients with hyperprogression using models 3 and 4 yielded AUCs of 0.877 and 0.774, respectively, with maximum gray value having great significance for the identification of hyperprogression.

At present, 3 methods are generally used to identify pseudoprogression and hyperprogression. The first method is biopsy of lesions. T cell recruitment and inflammatory cell infiltration have been confirmed as causes of pseudoprogression [8,9]. However, biopsy is an invasive procedure and it may be difficult to biopsy some deep lesions. The second method is response evaluation criteria for immunotherapy, including immune-related response criteria (irRC) [5]. and immune-related response evaluation criteria in solid tumors (irRECIST) [14]. These criteria allow PD patients who are generally eligible to continue immunotherapy and re-evaluate their responses at prescribed times (4 weeks or 4–8 weeks). The third method is tumor dynamic growth index. At present, however, definitions and diagnostic criteria have not been standardized. Indices used to evaluate tumor dynamic growth include tumor growth rate (TGR) [6], progress pace [7], and tumor growth kinetics (TGK) [15]. Each of these methods has limitations, making it impossible to accurately identify patients with pseudoprogression and hyperprogression. Using several filtered CT textures, entropy and consistency were found to be reduced after neoadjuvant chemotherapy for esophageal cancer, with the change in skewness being related to survival [16]. Moreover, 24 textures were found to predict lymph node metastasis in patients with colon cancer [17], and voxel heterogeneity analysis was better able to predict the curative effect of treatment in patients with rectal cancer than the traditional average volume analysis [18].

Radiomics has also shown promise in tumor diagnosis, clinical staging, response evaluation and prognosis prediction in patients with malignant tumors of the digestive system, including in patients receiving immunotherapy. For example, radiomics was able to predict response to immunotherapy by assessing tumor-infiltrating CD8 cells [19]. Moreover, 1860 radiomics features extracted at baseline could be used in patient classification and to predict those who will subsequently develop immunotherapy-induced pneumonitis [20]. Another study demonstrated that a CT-determined radiomics signature could identify tumors with increased lymphocyte infiltration and distinguish between pseudoprogression and true progression [21].

Radiomics has important advantages, including use of results of non-invasive imaging examinations. This allows tumors to be comprehensively and repeatedly evaluated. Imaging examinations provide not only general anatomical information but also functional information about the organ. Second, imaging results are traditionally interpreted by individual physicians, relying mainly on their subjective judgment. Thus, predictive ability is related to physician experience. In contrast, radiomics evaluates imaging results objectively and quantitatively, without the need for subjective determinations. Finally, research on the molecular mechanism of tumors has led to new treatment methods, such as targeted therapy and immunotherapy, which have greatly changed the prognosis of cancer patients. However, unlike traditional cytotoxic drugs, targeted therapy and immunotherapy may not induce tumor cell necrosis. Therefore, evaluations based on traditional measurements may not accurately reflect the anti-tumor effects of these agents. Appropriate screening of massive amounts of information is a key step in radiomics. Using the selected features to construct a model may allow multidimensional evaluation of tumor tissue.

The diagnostic ability of radiomics in patients with colorectal cancer and lung cancer may be improved by the inclusion of traditional clinical indicators, such as age, sex, and laboratory results, in the constructed model [22]. In practice, the latter factors, including demographic characteristics and the results of laboratory and pathological examinations, have great value in determining disease stage and predicting response to treatment. Therefore, the integration of multidimensional information into radiomics models could further improve their accuracy and efficiency.

Because radiomics is at the intersection between medicine and engineering, it requires collaborations among clinicians, imaging physicians, and data analysis professionals. Thus, use of radiomics in patient evaluation will require the participation of experts in these fields. Finally, radiomics is currently in the preclinical research stage, with most studies to date being single-center retrospective evaluations. Use of radiomics in clinical practice will require confirmation of the accuracy of radiomics models in large, multi-center studies.

Conclusion

This study explored the potential of radiomics in evaluating patient response to immunotherapy. The model constructed from CT features at baseline and first evaluation had the best classification performance, with a sensitivity, specificity, and AUC of 83.3%, 88.9%, and 0.806, respectively. The findings of this study suggest that radiomics can predict patient response to immunotherapy and may supplement conventional response evaluation criteria.