RNA sequencing profiles and diagnostic signatures linked with response to ramucirumab in gastric cancer.

Gastric cancer (GC) is the fifth-ranked cancer type by associated mortality. The proportion of early diagnosis is low, and most patients are diagnosed at the advanced stages. First-line therapy standardly includes fluoropyrimidines and platinum compounds with trastuzumab for HER2-positive cases. For recurrent disease, there are several alternative options including ramucirumab, a monoclonal therapeutic antibody that inhibits VEGF-mediated tumor angiogenesis by binding with VEGFR2, alone or in combination with other cancer drugs. However, overall response rate following ramucirumab or its combinations is 30%-80% of the patients, suggesting that personalization of drug prescription is needed to increase efficacy of treatment. We report here original tumor RNA sequencing profiles for 15 advanced GC patients linked with data on clinical response to ramucirumab or its combinations. Three genes showed differential expression in the tumors for responders versus nonresponders: CHRM3, LRFN1, and TEX15 Of them, CHRM3 was up-regulated in the responders. Using the bioinformatic platform Oncobox we simulated ramucirumab efficiency and compared output model results with actual tumor response data. An agreement was observed between predicted and real clinical outcomes (AUC ≥ 0.7). These results suggest that RNA sequencing may be used to personalize the prescription of ramucirumab for GC and indicate potential molecular mechanisms underlying ramucirumab resistance. The RNA sequencing profiles obtained here are fully compatible with the previously published Oncobox Atlas of Normal Tissue Expression (ANTE) data.

INTRODUCTION

Approximately 1 million patients are annually diagnosed with gastric cancer (GC), and nearly 740,000 die of this disease (Ferlay et al. 2010). Despite recent advances in GC treatment, long-term survival remains poor. Standard first-line therapy for unresectable, advanced, or metastatic gastric adenocarcinoma usually includes fluoropyrimidines and platinum compounds (Al-Batran et al. 2008), whereas more options are available at recurrence: ramucirumab (Fuchs et al. 2014), paclitaxel (Ajani et al. 1994), their combination (Wilke et al. 2014), docetaxel (Albertsson et al. 2007), and irinotecan (Thuss-Patience et al. 2011). Moreover, anti-PD1 immunotherapeutics can be considered an emerging treatment modality for metastatic GC (Brar and Shah 2019). In addition, molecular testing is required to define HER2 status of individual GC cases that guides prescription of trastuzumab (Bang et al. 2010). However, selection of other drugs including ramucirumab is not based on any molecular diagnostic features.

Ramucirumab is a monoclonal therapeutic antibody that inhibits VEGF-mediated tumor angiogenesis by binding with VEGFR2 (Krupitskaya and Wakelee 2009). For the recurrent/advanced/metastatic GC it may be used either as monotherapy (Smit et al. 2018; Yamaguchi et al. 2018) or in combinations with other drugs (e.g., paclitaxel) (Wilke et al. 2014) and therapeutic regimens like FOLFIRI (folinic acid + irinotecan + 5-fluorouracil) or FOLFOX (oxaliplatin + leucovorin + 5-fluorouracil) (Yoon et al. 2016) or cisplatin + capecitabine (Fuchs et al. 2019). In published clinical trials, overall response rate following treatment with ramucirumab or its above combinations was ∼30%–80% of the patients (Wilke et al. 2014; Yoon et al. 2016; Yamaguchi et al. 2018; Fuchs et al. 2019). The overall response rate includes a proportion of the patients who showed complete or partial response to the treatment or disease stabilization. Moreover, for the outcomes of combinational therapies, the impact of ramucirumab itself was unclear, thus suggesting that personalization of drug prescription is needed to increase efficacy of treatment (Yamaguchi et al. 2018).

Gene expression profiles can be considered attractive next-generation biomarkers for cancer molecular diagnostics (Buzdin et al. 2019b). Expressions of single genes (Patel and Kurzrock 2015) and their statistical combinations (van de Vijver et al. 2002) or activation levels of the whole molecular pathways (Buzdin et al. 2018) may be useful to predict individual tumor response to the therapies. In turn, RNA sequencing is currently accepted as the gold standard analytic approach for high-throughput screening of gene expression (SEQC/MAQC-III Consortium 2014). Comparison of tumor gene expression profiles with the normal tissues is crucial for understanding individual mechanisms of cancer development, progression, and response to the targeted therapies (Buzdin et al. 2019c; Rodon et al. 2019).

However, RNA sequencing data sets obtained using different equipment, reagents, and protocols may be poorly compatible with each other (Buzdin et al. 2014; Borisov et al. 2019), and ideally the same experimental platform should be used to compare the results (Borisov et al. 2019). We recently published an annotated database of RNA sequencing profiles termed Oncobox Atlas of Normal Tissue Expression (ANTE) (Suntsova et al. 2019) that represents 142 solid tissue samples from human healthy donors killed in road accidents and 17 blood samples from healthy volunteers. It has statistically significant reference groups for 20 human tissues/organs including the stomach.

In this study, we report clinically annotated RNA sequencing profiles for 15 advanced GC patients that are fully compatible with the ANTE collection of normal tissues because they were profiled using the same equipment, reagents, and protocols with the same analytic thresholds. Among the others, the molecular data are linked with the information on clinical response to ramucirumab or its combinations. We used bioinformatic platform Oncobox (Sorokin et al. 2018; Buzdin et al. 2019a; Poddubskaya et al. 2019a) to model ramucirumab efficiency based on RNA sequencing data and compared the output results with the known tumor response information. The Oncobox algorithm builds a personalized rating of target drugs for individual cancer patients. It is based on a simultaneous analysis of gene expression and molecular pathway activation in the patient's tumor and was shown to be effective in a prospective cohort of advanced cancer cases (Poddubskaya et al. 2018, 2019b) and in an ongoing prospective clinical investigation (Poddubskaya et al. 2019a). We observed here a statistically significant coincidence of the predicted and real GC clinical outcomes (AUC ≥ 0.7). These results demonstrate that gene expression analyses may uncover molecular mechanisms underlying ramucirumab resistance and, more importantly, may be useful to personalize prescription of ramucirumab to recurrent GC patients.

RESULTS

In this study, we investigated links between experimental gene expression profiles established by RNA sequencing of GC biosamples and the patient's clinical responses to the therapy containing therapeutic monoclonal antibody ramucirumab.

Clinical Data

In this study, we enrolled 15 patients with advanced or metastatic GC. The patients were receiving ramucirumab-containing therapies in two oncology hospitals: M.F. Vladimirsky Moscow Regional Research Clinical Institute and N.N. Blokhin National Medical Research Center of Oncology. In total, there were five female and 10 male patients enrolled with mean ages of 47.2 and 57.6 yr old, respectively (Table 1). All patients had progressive disease before receiving ramucirumab-containing therapy. Ramucirumab was prescribed as second-line therapy for both metastatic and primary tumors according to the National Comprehensive Cancer Network (NCCN) guidelines for stomach cancer (https://www.nccn.org/patients/guidelines/content/PDF/stomach.pdf). In 14 cases, primary GC specimens were analyzed, whereas in one case, a peritoneal metastasis of GC was investigated. The GC histological subtypes were tubular adenocarcinoma (n = 8), poorly cohesive adenocarcinoma, signet-ring cell type (n = 5), poorly cohesive adenocarcinoma, not otherwise specified (NOS) (n = 1), and omental metastasis of mucinous gastric adenocarcinoma (n = 1) (Table 2; Fig. 1). PCA analysis indicated that normal stomach and brain tissues from the ANTE database (Suntsova et al. 2019) formed tight clusters, which were not mixed with GC samples from this study (Fig. 1E). In addition, GC samples were closer to normal stomach than to normal brain, indicating compatibility of the data sets (Fig. 1E). We used normal brain tissues as an outgroup in order to test the hypothesis that cancer and normal data sets in this study are compatible. Compatibility of data sets implies that biological differences (that manifest themselves as variance in gene expression space) are higher than any of the between data set differences. A PCA plot (Fig. 1E) shows that biological differences between gastric and neural tissues are higher than any of differences between two data sets. All patients underwent ramucirumab therapy either as monotherapy (n = 7) or in combination with paclitaxel (n = 6) or FOLFIRI regimen (n = 2) (Table 3). The registered clinical outcomes of treatment were “partial response,” “stable disease,” and “progressive disease.” In this study, the patients were classified as either responders—for “partial response” and “stable disease” outcomes—or nonresponders—for “progressive disease” outcomes (Table 3).

Table 1.

Patients’ demographic and clinical information

Patient ID	Sex	Age	Primary tumor/metastasis	TNM	Number of lines of therapy prior to ramucirumab
GC_4	Male	67	Primary	T4aN3M0	1
GC_5	Male	45	Primary	T4N1M1	4
GC_6	Male	65	Primary	T4bNxM1	1
GC_7	Male	79	Primary	T3N0M0	2
GC_8	Female	51	n/a	T3NxM1	2
GC_9	Male	59	Primary	T4aN3bM0	2
GC_11	Female	40	Primary	TxNxM1	1
GC_12	Male	41	Primary	T3NxM0	1
GC_13	Female	54	Primary	T4aNxM1	1
GC_14	Male	33	Primary	T4N1M1	3
GC_15	Female	31	Primary	T4N2M1	1
GC_16	Male	57	Primary	T4N2M0	1
GC_17	Female	60	Peritoneum metastasis	T4bN2M1	1
GC_18	Male	63	Primary	T4bN3M1	1
GC_19	Male	67	Primary	T4N1M1	1

Table 2.

Histological typing of gastric tumor samples investigated

Patient ID	Histology
GC_4	Poorly cohesive adenocarcinoma, signet-ring cell type
GC_5	Tubular adenocarcinoma
GC_6	Poorly differentiated tubular adenocarcinoma
GC_7	Tubular adenocarcinoma
GC_8	Poorly cohesive adenocarcinoma, signet-ring cell type
GC_9	Poorly differentiated tubular adenocarcinoma
GC_11	Poorly cohesive adenocarcinoma, signet-ring cell type
GC_12	Poorly differentiated tubular adenocarcinoma
GC_13	Poorly cohesive adenocarcinoma, signet-ring cell type
GC_14	Tubular adenocarcinoma
GC_15	Poorly cohesive adenocarcinoma, signet-ring cell type
GC_16	Poorly cohesive adenocarcinoma NOS
GC_17	Omental metastasis of mucinous gastric adenocarcinoma
GC_18	Tubular adenocarcinoma
GC_19	Poorly differentiated tubular adenocarcinoma

Figure 1.

Histological subtypes of gastric cancer (GC) samples and principal component analysis (PCA). (A) Tubular adenocarcinoma; (B) poorly cohesive adenocarcinoma, signet-ring cell type; (C) poorly cohesive adenocarcinoma, NOS; (D) omental metastasis of mucinous gastric adenocarcinoma; and (E) PCA of stomach cancer samples, normal stomach samples from the ANTE database, and normal brain samples (outgroup) from the ANTE database (Suntsova et al. 2019).

Table 3.

Ramucirumab treatment regimen for patients enrolled and their clinical outcome

Patient ID	Treatment regimen	Clinical outcome	PFS (months)
GC_4	Paclitaxel (80 mg/m2 days 1, 8, 15) + ramucirumab (8 mg/kg days 1, 15) − 28 d cycles, followed by 12 cycles of ramucirumab monotherapy (ramucirumab 8 mg/kg days 1, 15) − 28 d cycles	Stable disease	11
GC_5	Ramucirumab monotherapy	Progressive disease	1
GC_6	Seven courses: FOLFIRI (fluorouracil + leucovorin + irinotecan) + ramucirumab (8 mg/kg days 1, 15)	Stable disease	5
GC_7	Ramucirumab monotherapy	Progressive disease	5
GC_8	Two courses paclitaxel + ramucirumab	Progressive disease	2
GC_9	Paclitaxel + ramucirumab	Progressive disease	9
GC_11	Two courses: paclitaxel (80 mg/m2 days 1, 8, 15) + ramucirumab (8 mg/kg days 1, 15)	Progressive disease	3
GC_12	Ramucirumab (8 mg/kg every 14 d)	Progressive disease	6
GC_13	Eight FOLFIRI + ramucirumab (8 mg/kg every 14 d)	Stable disease	8
GC_14	Paclitaxel + ramucirumab	Progressive disease	1
GC_15	11 injections: ramucirumab 8 mg/kg	Partial response	8
GC_16	15 injections: ramucirumab 8 mg/kg	Stable disease	7
GC_17	17 injections: ramucirumab 8 mg/kg + 4 injections paclitaxel	Stable disease	11
GC_18	14 injections: ramucirumab 8 mg/kg	Stable disease	8
GC_19	Six injections: ramucirumab 8 mg/kg	Progressive disease	5

RNA Sequencing and Molecular Data Quality Control

The GC tissue specimens were the formalin-fixed paraffin-embedded (FFPE) tissue blocks stored in a clinical diagnostic laboratory for 4–49 mo before extraction of RNA. Previously we used the same experimental protocol for generating the Oncobox atlas of RNA sequencing profiles for normal human tissues from healthy donors (Suntsova et al. 2019). We found that the profiles with the number of sequencing reads mapped on human genes exceeding 2.5 million were clustering tissue-specifically, whereas for the profiles with a lower number of mapped reads we observed biased clustering. A quality control (QC) threshold of 2.5 million gene-mapped reads was, therefore, established to identify high-quality RNA sequencing profiles (Suntsova et al. 2019). In this study, we used the same Oncobox RNA sequencing and data processing protocol. All the current 15 GC RNA sequencing profiles passed the QC threshold and were further investigated. The original sequencing data were deposited to the National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA) with accession number PRJNA562149. Summary statistics of RNA sequencing (number of unmapped, uniquely mapped, and multimapped reads as well as average exon coverage) are depicted in Supplemental Table S1.

Differentially Regulated Genes and Molecular Pathways

According to differentially expressed gene analysis using DESeq2 software, only three genes passed the adjusted P-value threshold of 0.05 for “responders” on ramucirumab-containing therapies versus “nonresponders.” These genes were CHRM3 (cholinergic receptor muscarinic 3), LRFN1 (leucine rich repeat and fibronectin type III domain containing 1), and TEX15 (testis expressed 15, meiosis and synapsis associated). Of them, CHRM3 was up-regulated in the responders, whereas the other genes were down-regulated. We found no previous literature reports on the implication of these genes in GC. However, CHRM3 is involved in up-regulation of MAPK signaling and invasion and migration of prostate and colorectal cancer cells (Belo et al. 2011; Zheng et al. 2019), and its genetic variants are associated with polycystic ovary syndrome (Kim et al. 2019) and bladder cancer (Wang et al. 2016) in Korean and Chinese populations. It is also known as the poor-prognosis biomarker in endometrial carcinoma (Wang et al. 2015). In turn, mutations in the DNA repair gene TEX15 are known to be linked with a high risk of prostate and breast cancers (Lin et al. 2017).

We then calculated pathway activation levels of 3125 molecular pathways using the Oncobox software (Sorokin et al. 2018) and tested them for differential activation between the responder and nonresponder tumors (Supplemental Table S2). None of the pathways passed FDR threshold of 0.05, but the most significantly differential pathway according to the Wilcoxon rank sum test was the “Nectin adhesion pathway (positive regulation of JNK cascade).” This pathway is a fragment of the Nectin adhesion regulatory network, which is responsible for downstream positive regulation of JNK (c-Jun amino-terminal kinase) cascade. This pathway appeared to be less active in the treatment responders compared with the nonresponders, mostly as a result of decreased expression of RAP1A, RAP1B, and SRC, which are upstream regulators of the CDC42 gene product (Fig. 2). Interestingly, activation of JNK signaling in GC cells can lead to enhanced resistance against platinum-based chemotherapeutics (Ye et al. 2015) and microtubule-targeting drugs including paclitaxel (Cui et al. 2017).

Figure 2.

The nectin adhesion pathway (positive regulation of JNK cascade) was more active in the nonresponder tumors (upper panel) than in the responder tumors (lower panel). The pathway is shown as an interacting network in which green arrows indicate activation. Gene expression values were geometrically averaged for responders and nonresponders. The depth of color of each node of the network corresponds to the logarithms of the case-to-normal (CNR) expression rate for each node, in which “normal” is a geometric average between normal tissue samples. The normal tissue gene expression profiles for healthy donors were taken from the ANTE database (Suntsova et al. 2019). The scale represents the extent of up/down-regulation of the pathway nodes. PAL (pathway activation level) was calculated using the Oncobox software (described in Methods).

To assess the variance of the pathway activation level (PAL) and gene expression in the Nectin adhesion pathway (positive regulation of JNK cascade) between the cohorts of responders and nonresponders, we performed a nonpaired t-test and Wilcoxon test. We obtained for the t-test P = 0.01, for the Wilcoxon test P = 0.0003, and an absolute difference of average PAL scores of 0.35, thus indicating that the observed difference is significant considering variances of PAL in each cohort. Similarly, we performed DESeq2 differential gene expression analysis for all genes forming this pathway (Table 4). For all individual genes, the expression differences did not reach a P < 0.05 significance threshold, thus confirming previous findings on superior robustness of molecular pathways as cancer biomarkers (Borisov et al. 2014, 2017; Ozerov et al. 2016; Buzdin et al. 2018).

Table 4.

DESeq2 differential analysis results for the five genes of the “Nectin adhesion pathway (positive regulation of JNK cascade)” pathway

Gene	Log2 (fold change)	P-value
RAP1A	−0.81	0.12
CDC42	−0.67	0.17
RAP1B	−0.34	0.21
SRC	−0.35	0.22
FARP2	−0.05	0.84

To our knowledge, none of those three differential genes (TEX15, CHRM3, LRFN1)/molecular pathways identified was previously connected with the sensitivity to ramucirumab in cancer cells.

We then used the Oncobox platform to build a personalized rating of targeted therapeutics for each individual GC sample according to their simulated abilities to inhibit aberrantly regulated molecular pathways and drug target genes. To this end, the RNA sequencing profiles obtained for the GC samples were individually compared with the set of profiles for normal stomach tissues obtained from the ANTE database (Suntsova et al. 2019). The output modeled drug efficiency values, termed balanced efficiency scores (BESs), which were obtained for 159 targeted therapeutics and prioritized drugs according to this score. A higher BES suggests higher predicted efficacy of a drug for an individual tumor (Poddubskaya et al. 2019b). We compared Oncobox BES values with the actual clinical outcomes of 15 ramucirumab-treated patients (Supplemental Table S3) and found that BES values could predict the GC treatment response on ramucirumab-containing therapies with the quality metric area under the receiver operating characteristic (ROC) curve (ROC AUC) of 0.7 (Fig. 3A).

Figure 3.

The receiver operating characteristic (ROC) area under the curve (AUC) and histograms for predicting response to ramucirumab in gastric cancer (GC) patients using Oncobox balanced efficiency score (BES) values. (A) ROC AUC for the complete set of patients (n = 15). (B) ROC AUC for the patients who received ramucirumab monotherapy (n = 7). (C) Histogram and threshold value for the complete set of patients (n = 15). (D) Histogram and threshold value for the patients who received ramucirumab monotherapy (n = 7). The Oncobox BES of ramucirumab was used as a predictor of response. Patients were divided into two classes: “nonresponders,” for progressive disease, and “responders,” for stable disease and partial response (RECIST). Threshold values were assigned based on an assumption of the equality of type I and type II error rates.

The AUC value is the universal characteristic of biomarker robustness depending on its sensitivity and specificity (Green and Swets 1966). This statistical approach is broadly applicable to different types of biomarkers in oncology (Borisov et al. 2014; Lezhnina et al. 2014; Liu et al. 2018; Tanioka et al. 2018; Zolotovskaia et al. 2018; Chen et al. 2019). It positively correlates with the quality of a biomarker and varies from 0.5 to 1. The standard discrimination threshold is 0.7, and the parameters with greater AUC are considered good-quality biomarkers, and vice versa (Boyd 1997).

The AUC value of 0.7 obtained for the whole set of patients (Fig. 3A) therefore suggests that the Oncobox BES value may serve as the good-quality predictive biomarker of GC clinical response for patients on ramucirumab. Importantly, reducing the patient cohort to those who received only ramucirumab monotherapy (seven cases) significantly improved performance of the BES-based predictions by increasing AUC to 0.75 (Fig. 3B). This can be explained by interference of other drugs (FOLFIRI, paclitaxel) in the non-monotherapy regimens.

To further study statistical performance of BES values, we pick up threshold values for cohorts of all patients (15 cases) and of those who received only ramucirumab monotherapy (seven cases). These thresholds were introduced to separate patients with different response status assuming equal importance for type I and type II errors (Fig. 3C,D); confusion matrices are shown in Tables 5 and 6. The sensitivity and specificity values observed were, respectively, 0.83 and 0.71 for all patients, and 1 and 0.75 for the monotherapy cohort.

Table 5.

Confusion matrix for a cohort of all patients (15 cases) using the threshold BES = 0.1, with sensitivity of 0.62 and specificity of 0.14

	Cases predicted as responders	Cases predicted as nonresponders
Responders	5	3
Nonresponders	6	1

Table 6.

Confusion matrix for a cohort of all patients (15 cases), with sensitivity of 0.83 and specificity of 0.71

	Cases predicted as responders	Cases predicted as nonresponders
Responders	6	1
Nonresponders	2	5

In an ongoing prospective trial of the Oncobox platform as a second-opinion tool for prescription of targeted therapeutics (NCT03724097), a BES threshold of 0.1 was arbitrarily set to discriminate between potentially effective and ineffective drugs. In this study, we tested this threshold for the original patient survival data. Eleven patients had greater and four patients had lower BES values (Fig. 4). Progression-free survival analysis revealed that Oncobox BES was able to predict better treatment outcomes for ramucirumab-containing therapies with hazard ratio 0.16 (95% CI, 0.034–0.74, P = 0.019; Fig. 4) and log-rank test P = 0.013. However, Cox proportional hazards analysis returned Cox P = 0.31. Such discrepant results between log-rank and Cox proportional hazards tests can be caused by small sample size. However, our further analysis showed that this BES threshold of 0.1 was suboptimal in terms of sensitivity and specificity values obtained, which were 0.62 and 0.14, respectively (confusion matrix shown on Table 5). Our results, therefore, suggest that using a greater BES threshold can be beneficial for algorithmic discrimination of the potential treatment responders and nonresponders, at least for the case of using ramucirumab in GC. Here we exemplify that a threshold of 3.1 results in much better values of sensitivity and specificity for our group of samples (0.83 and 0.71, respectively) (Table 6; Fig. 3C). Taken together, these data strongly suggest that the individual tumor responses for patients on ramucirumab can be predicted by using RNA sequencing profiles as the complex set of biomarkers.

Figure 4.

The Oncobox BES for ramucirumab predictive power for progression-free survival in GC patients receiving chemotherapy regimens with ramucirumab (N = 15). Cohorts of patients with BES values lower or higher than 0.1 are shown separately. The risk table is shown below the plot. Hazard ratio = 0.16 (95% CI, 0.034–0.74, P = 0.019). Log-rank test P = 0.013.

DISCUSSION

Ramucirumab is the therapeutic monoclonal anti-VEGFR2 antibody that is approved for advanced GC or gastroesophageal junction adenocarcinoma either as monotherapy or in combination with paclitaxel (Fuchs et al. 2014; Wilke et al. 2014). It has been also approved for hepatocellular carcinoma (Zhu et al. 2019), non-small-cell lung cancer (Garon et al. 2014), and colorectal cancer (Tabernero et al. 2015). However, to our knowledge, so far there are no accepted biomarkers of ramucirumab efficiency in cancers (Bignucolo et al. 2017; Vlachostergios et al. 2018). We were also unable to find in the public domain any clinically annotated high-throughput gene expression profiles linked with the patient responses on treatment with ramucirumab. Here we for the first time experimentally analyzed RNA sequencing profiles of GC tumors with the available clinical outcomes of ramucirumab-containing treatment regimens.

When comparing treatment responders versus nonresponders, we found three statistically significantly differentially expressed genes: CHRM3, LRFN1, and TEX15. The most significantly differential molecular pathway was the “Nectin adhesion pathway (positive regulation of JNK cascade).” It is a fragment of the Nectin adhesion regulatory network that is responsible for positive regulation of the JNK cascade. We found that this pathway was more active in the nonresponders (Fig. 2) as a result of stronger expression of the RAP1A, RAP1B, and SRC genes, which encode upstream regulators of the CDC42 protein. CDC42 is a member of the Rho GTPase family that is involved in regulation of cell migration in various types of human cancers (Maldonado and Dharmawardhane 2018). In particular, CDC42 inhibition significantly decreased proliferation, migration, and invasion of GC cells (Du et al. 2016). Therefore, down-regulation of CDC42 activators, such as RAP1A, RAP1B, and SRC, may potentially be linked with tumor response to ramucirumab, because ramucirumab is an anti-VEGFR2 antibody that inhibits various signaling pathways including the SRC pathway (Clarke and Hurwitz 2013). Thus, decreased activation of the Nectin adhesion pathway (positive regulation of JNK cascade) in tumors may be linked to a response to ramucirumab in GC tissues.

Targeting molecular pathways abnormally activated in individual tumors was previously proposed to be a promising strategy for personalized selection of targeted therapies (Buzdin et al. 2018). Applicability of this approach to real-world cancer cases was tested in several clinical reports (Poddubskaya et al. 2018, 2019b). Here we present a retrospective cohort of advanced GC patients with known clinical outcomes for ramucirumab-containing therapies linked with high-throughput RNA sequencing profiles. Gene expression data were analyzed using the Oncobox bioinformatical platform for modeling the patient's tumor tissue sensitivity to ramucirumab. The predicted efficiencies were correlated with the actual tumor responses and progression-free survival data. Our results suggest that success of targeted therapies could be enhanced if guided by personalized analysis of gene expression. In this study, we experimentally profiled and bioinformatically investigated molecular data for 15 GC patients, with the smallest predicted response group including only four patients (ramucirumab BES < 0.1). As far as we know, this is currently the largest group of publicly available high-throughput gene expression profiles associated with known response statuses of patients being treated with ramucirumab. However, larger patient cohorts will be required to validate the trends observed here and to develop transcriptomics-based companion diagnostics for ramucirumab to personalize its use in GC.

METHODS

Tissue Samples

All experimental biosamples were FFPE tumor tissue blocks. All biosamples were evaluated by a pathologist to confirm the tumor tissue origin, and only the specimens with tumor cell content of >50% were further investigated. Of them, five GC samples were obtained from M.F. Vladimirsky Moscow Regional Research Clinical Institute and 10 GC samples were obtained from N.N. Blokhin National Medical Research Center of Oncology. The patients were 10 men and five women, whose mean ages were 57.6 and 47.2 yr old, respectively (with a range of 31–79 yr old).

The samples were clinically annotated and contained information about the patient's sex, age, diagnosis, and clinical history. For all the biosamples, informed written consent to participate in this study was collected from the patients or their legal representatives. The consent procedure and the design of the study were approved by the ethical committee of Vitamed Clinic, Moscow, protocol date 16.10.17. Clinical annotation of the tumor tissues investigated is summarized in Tables 1–3.

Preparation of Libraries and RNA Sequencing

To isolate RNA, 10-µM-thick paraffin slices were trimmed from each FFPE tissue block using a microtome. RNA was extracted from FFPE slices using QIAGEN RNeasy FFPE Kit following the manufacturer's protocol. RNA 6000 Nano or Qubit RNA Assay kits were used to measure RNA concentration. The RNA integrity number (RIN) was measured using the Agilent 2100 Bio-Analyzer. For depletion of ribosomal RNA and library construction, the KAPA RNA Hyper with rRNA erase kit (HMR only) was used. Different adaptors were used for multiplexing samples in one sequencing run. Library concentrations and quality were measured using the Qubit ds DNA HS Assay kit (Life Technologies) and Agilent TapeStation. RNA sequencing was done at the Department of Pathology and Laboratory Medicine, University of California Los Angeles, using Illumina HiSeq 3000 equipment for single-end sequencing, 50-bp read length, for ∼30 million (mln) raw reads per sample. A data quality check was done on Illumina SAV. Demultiplexing was performed with the Illumina Bcl2fastq2 v 2.17 program. Sequencing data were deposited in the NCBI SRA under accession ID PRJNA562149. Summary statistics of RNA sequencing (number of unmapped, uniquely mapped, and multimapped reads as well as average exon coverage) is depicted in Supplemental Table S1.

Processing of RNA Sequencing Data

RNA sequencing FASTQ files were processed with STAR aligner (Dobin et al. 2013) in “GeneCounts” mode with the Ensembl human transcriptome annotation (Build version GRCh38 and transcript annotation GRCh38.89). Ensembl gene IDs were converted to HGNC gene symbols using the Complete HGNC data set (https://www.genenames.org/; database version from 2017 July 13). Expression levels were established for 36,596 annotated genes with the corresponding HGNC identifiers. The minimum number of uniquely mapped reads was 3.75 mln for studied biosamples, with the mean value of 9.78 mln. Differential gene expression analysis was performed using DESeq2 software (Love et al. 2014). PCA was performed for log-transformed DESeq2 normalized counts using the R prcomp function.

Molecular Pathway Analysis and Ranking of Target Drugs

Pathway activation levels were established using the Oncobox analytic software (Sorokin et al. 2018) for 3125 molecular pathways extracted from the public databases Reactome (Croft et al. 2014), NCI Pathway Interaction Database (Schaefer et al. 2009), Kyoto Encyclopedia of Genes and Genomes (Kanehisa and Goto 2000), HumanCyc (Romero et al. 2005), Biocarta (Nishimura 2001), and QIAGEN pathway-central (available at https://www.qiagen.com/us/shop/genes-and-pathways/pathway-central/). The pathways investigated here are listed in Supplemental Table S1. The molecular pathways were visualized using the Oncobox pathway visualization/reconstruction tool (Sorokin et al. 2018; Buzdin et al. 2019a).

The PAL scores were calculated according to the Oncobox method (Borisov et al. 2020). For the pathways interrogated, gene expression values were geometrically averaged for the treatment responder and separately for the nonresponder biosamples. These averaged gene expression data were normalized on the normal tissue expression profiles extracted from the ANTE database (Suntsova et al. 2019). PAL scores were calculated as follows: in which PAL is pathway activation level, n is number of genes in a pathway, average_gene_expression_cancer is geometrically averaged gene expression value in all samples of a given cohort (treatment responders or nonresponders), and average_gene_expression_norm is geometrically averaged normalized gene expression value in control samples. Ranking of target cancer drugs using BES was performed as described previously (Buzdin et al. 2018; Poddubskaya et al. 2019b; Tkachev et al. 2020). The Oncobox software returned a personalized list of target drugs in descending order of predicted efficacy. The observed clinical responses were used for validation of the Oncobox predictions using ROC AUC analysis. ROC AUC was calculated using the R ROCR package. Patient survival analysis was performed using the R ggsurvplot package.

ADDITIONAL INFORMATION

Data Deposition and Access

The original sequencing data were deposited to NCBI SRA (https://www.ncbi.nlm.nih.gov/sra) with accession number PRJNA562149, allowing free access.

Ethics Statement

For all the biosamples, informed written consent to participate in this study was collected from the patients or their legal representatives. The consent procedure and the design of the study were approved by the ethical committee of Vitamed Clinic, Moscow, protocol date 16.10.17.

Acknowledgments

We acknowledge the patients for their participation in this study.

Author Contributions

M.Se., E.P., M.So., A.Ga., and A.B. conceived of and designed the study. A.Gl., E.P., M.B., L.K., and E.M. collected and interpreted patient data. M.Se., M.Su., L.K., E.M., M.B., and E.P. isolated and prepared tissue samples. D.A. did histological analysis. M.Su. and D.A. performed molecular analyses. M.So., V.T., and A.B. analyzed the data. M.Se., A.B., and M.So. wrote the manuscript.

Funding

This study was supported by the OmicsWay research program in oncology, and by Russian Science Foundation grant 18-15-00061 (M.So., M.Su., and A.B.).

Competing Interest Statement

M.So., A.Ga., V.T., and A.B. have a financial relationship with Omicsway Corp. The remaining authors have nothing to disclose.