Predictive biomarkers of immunotherapy response with pharmacological applications in solid tumors

1Immune-checkpoint inhibitors (ICIs) were introduced for the treatment of solid and hematological malignancies with outstanding results in the last decade. There are three groups of ICIs. The first group consists of ICIs inhibiting cytotoxic T lymphocyte-associated protein 4 (CTLA-4) on T cells, the second one is related to programmed cell death 1 (PD-1) receptor on lymphocytes, and the third group is linked to programmed cell death ligand 1 (PD-L1) on tumor cells. Physiologically, PD-1 is expressed on several immune cells (e.g., lymphocytes, natural killer cells), and PD-L1 is present on almost all somatic cells (e.g., hematopoietic cells or tumor cells). In tumors, blockade of the PD-1/PD-L1 axis results in the survival of the malignant cells.

2The first ICI approved by the U.S. Food and Drug Administration (FDA) for treating metastatic melanoma was ipilimumab, a fully human monoclonal antibody against CTLA-4. Ipilimumab was shown to be effective in an extremely wide range of advanced cancers including melanoma, renal cell carcinoma (RCC), mismatch repair deficient (dMMR) or microsatellite instability-high (MSI-H) colorectal carcinoma, hepatocellular carcinoma (HCC), pleural mesothelioma, and non-small-cell lung cancer (NSCLC). Tremelimumab, another fully human monoclonal antibody against CTLA-4, has been investigated in multiple solid tumors with mixed results. Tremelimumab has not yet received an FDA-approval so far.

3Nivolumab was the first approved monoclonal antibody targeting PD-1 and can be administered in advanced melanoma, Hodgkin’s lymphoma, HCC, NSCLC and SCLC, RCC, head and neck cancer, urothelial carcinoma, CRC, gastric, or esophageal adenocarcinoma, and malignant pleural mesothelioma. Pembrolizumab, a humanized monoclonal antibody against PD-1 has been approved either as a monotherapy or as a combination therapy for the treatment of recurrent or metastatic melanoma, NSCLC, head and neck squamous cell cancer (HNSCC), gastric/gastroesophageal junction adenocarcinoma, lymphoma, urothelial cancer, cervical cancer, Merkel cell carcinoma, RCC, triple-negative breast cancer (TNBC), cutaneous squamous cell carcinoma, endometrial cancer, and HCC. Pembrolizumab can be administered to patients regardless of their age and tumor type in case their tumor is MSI-H or dMMR. Recently, dostarlimab, a humanized anti-PD-1 monoclonal antibody, has also received an FDA-approval for recurrent or advanced MSI-H/dMMR endometrial cancer and other solid tumors.

4PD-L1-blocking antibodies started with avelumab in 2015 for metastatic Merkel cell carcinoma, and then it continued with locally advanced or metastatic urothelial carcinoma (mUC), and advanced RCC. Durvalumab, another human monoclonal antibody was accepted for metastatic urothelial bladder cancer, urothelial carcinoma, NSCLC, and SCLC. Atezolizumab, a third PD-L1 monoclonal antibody, is also accepted for locally advanced or mUC, metastatic NSCLC, SCLC, TNBC, HCC, and melanoma.

5Biomarkers are biological indicators that can be used to select patients for a systemic anticancer treatment like immunotherapy. A major limitation of the widespread use of immunotherapies is the fact that there are only a few clinically useful biomarkers capable to predict therapy response. A study found that tumor mutational burden (TMB), and PD-L1 expression can predict response to pembrolizumab in a huge variety of cancers (melanoma, bladder cancer, breast cancer, CRC, HNSCC, and SCLC). The correlation between PD-L1 expression and MSI and response to pembrolizumab was also investigated in gastric cancer. Cluster of Differentiation 8 positive (CD8+) T cell phenotype and TMB were associated with enhanced response to atezolizumab in mUC. Findings from other studies highlighted the central role of the tumor microenvironment (TME) in nivolumab, pembrolizumab, and anti-CTLA-4 response. The importance of both innate, and adaptive immune systems was emphasized in connection with anti-PD-1 response in NSCLC, melanoma, and HNSCC. Meanwhile, in recent years, the tumors of several ICI-treated patient cohorts were investigated with transcriptomic technologies. The simultaneous analysis of the entire transcriptome makes it possible to identify new genes capable to serve as biomarkers of response and to validate previously suggested biomarker candidates.

6Here, our goal was to expose and integrate available transcriptomic datasets of ICI-treated tumors to establish a framework enabling an integrated analysis of genes related to treatment sensitivity or resistance. Uncovering robust genes with increased expression in treatment-resistant tumors could offer the opportunity to develop a targeted therapy to augment the effects of immune-checkpoint inhibitors.

7We screened the NCBI Gene Expression Omnibus (GEO) repository using the keywords “human [organism] AND (anti-PD-1 OR anti-PD-1 OR anti-PD-L1 OR anti-PD-L1 OR anti-CTLA-4 OR anti-CTLA-4)”, and “human [organism] AND (pembrolizumab OR nivolumab OR atezolizumab OR durvalumab OR avelumab OR cemiplimab OR ipilimumab OR camrelizumab OR cintilimab OR tislelizumab OR toripalimab)” on 10th Jan 2022. We omitted datasets with no available gene expression or clinical data, or with single-cell RNA-sequencing (scRNA-Seq), T or B cell receptor sequencing (TCR/BCR-Seq), non-mRNA-sequencing (e.g., whole-exome sequencing, non-coding RNA profiling, methylation profiling, protein array), studies of cell lines, stem cells, sorted peripheral blood mononuclear cells, studies in mice, studies without cancer, and GEO SuperSeries files. We also conducted a literature research to find additional studies. Our scope of investigation only included studies with simultaneously available clinical (response) and bulk-tissue gene expression data.

8Gene expression data from all eligible datasets were combined into a single table, quantile normalized and scaled to 1000. For the clinical annotation, we categorized patients as responders or non-responders based on reported pathological response or survival time. Those patients were selected as responders who experienced progression-free survival (PFS) longer than 12 months or had a partial response (PR) or complete response (CR). Those who experienced less than 12 months of PFS or had progressive disease (PD) or stable disease (SD) were categorized as non-responders. Survival time was not used if the patient had no event and the follow-up time was censored before 12 months. Tumor samples obtained before induction of the therapy were termed “pre-treatment” samples, and tumors collected during or after the therapy were termed “on-treatment” samples.

9Three separate analyses were performed across all genes to identify pre-treatment gene expression changes related to response against anti-PD-1, anti-PD-L1, and anti-CTLA-4 treatment. On-treatment samples were left out of the analysis because of low sample sizes.

10We used Gene Ontology (GO) Enrichment analysis to uncover biological processes connected to the gene lists related to response against anti-PD-1, anti-PD-L1, and anti-CTLA-4 treatment.

11Druggability of candidate genes was determined by a literature search in PubMed and GeneCards (https://www.genecards.org/) to include those where (1) in silico prediction, (2) in vitro assay, (3) clinical study or (4) FDA-approval of the given drug were available.

12We extended our previously established ROC plotter platform to enable the investigation of new biomarkers and the validation of current results in all patients treated with immunotherapy. The platform is running on Ubuntu 20.04.4 LTS server driven by Apache 2.4.41. The front-end site was developed in PHP using the YII2 framework. The application data are stored in the PostgreSQL database and the computations are performed via an R script. The portal can be reached at https://www.rocplot.com/immune.

13Statistical analysis was conducted in the R environment (https://www.r-project.org/) using Bioconductor libraries (https://www.bioconductor.org/). To find differentially expressed genes associated with improved or worse outcomes, Mann-Whitney unpaired U-test and receiver operating characteristic (ROC) analysis were used. To avoid false discovery due to multiple testing, Bonferroni-adjustment (P = 0.05) was applied with the service https://www.multipletesting.com/.

14We have identified 225 series files in NCBI GEO fulfilling the initial search criteria. Through literature research, we also found the Cancer Research Institute iAtlas (CRI iAtlas) (https://www.cri-iatlas.org/), another portal with ICI-treated samples, in which another six datasets were found. Finally, five additional cohorts were found by looking up the referenced literature.

15In summary, from NCBI GEO, CRI iAtlas, Chen et al., Litchfield et al., Liu et al., Kim et al., and Miao et al., altogether 246 datasets with 3823 samples were found. Then, we removed duplicate records. For example, Litchfield et al. described eleven studies out of which four were also included in CRI iAtlas, and two in GEO (GSE78220, GSE91061). A detailed description of the complete screening process is provided in Fig. 1.

16Datasets were individually investigated to include only those in which at least one of the following clinical variables was reported: progression-free survival or interval (PFS/PFI), relapse-free survival (RFS), overall survival (OS), recurrence, and response determined by the Response Evaluation Criteria in Solid Tumors (RECIST) (including complete response (CR), partial response (PR), stable disease (SD), or progressive disease (PD)).

17Twenty datasets with 2222 samples met these eligibility criteria and were selected for further analysis (Fig. 1). Sixty-eight samples were excluded due to different restrictions, such as (1) a very low number of samples treated with a specific agent (e.g., avelumab or experimental drugs), (2) duplicated samples, (3) missing clinical variable (e.g., event for survival), (4) no expression data, (5) ambiguous treatment time, (6) samples taken from metastatic sites or (7) low incidence of the given tumor type. By using these filtering criteria, the cohort was reduced to 2,154 samples out of which 720 samples were blood samples – these were removed from further analysis.

18Finally, 1434 samples from 19 datasets were investigated manually. Different tumor types were integrated into the database, including (a) melanoma (n = 570), (b) esophageal and gastroesophageal junction adenocarcinoma (n = 103), (c) gastric cancer (n = 45), (d) urothelial cancer (bladder/ureter/pelvis cancer) (n = 449), (e) head and neck squamous cell carcinoma (n = 110), (f) hepatocellular carcinoma (n = 22), g) lung cancer (small cell and non-small-cell lung cancer, or squamous and non-squamous non-small cell lung cancer) (n = 60), (h) breast cancer (triple-negative (n = 12), and ER + HER2- breast cancer (n = 2)), (i) renal cell carcinoma (n = 33), and (j) glioblastoma (n = 28). Patients either received anti-PD-1 (nivolumab, or pembrolizumab) (all n = 877; pre-treatment n = 776; on-treatment n = 101), anti-PD-L1 (atezolizumab, or durvalumab) (all n = 488; pre-treatment n = 457; on-treatment n = 31), or anti-CTLA-4 (ipilimumab) (all n = 124; pre-treatment n = 98; on-treatment n = 26) treatments, out of which 55 patients received combination therapy with anti-PD-1 (Fig. 1).

19Eventually, datasets were stratified into six groups: (1) pre-treatment anti-PD-1 (GSE78220, GSE91061, GSE93157, GSE115821, GSE121810, GSE136961, GSE140901, GSE165745, GSE176307, Cristescu2018, Gide2019, Kim2018, Miao2018, and Liu2019), (2) on-treatment anti-PD-1 (GSE78220, GSE91061, GSE115821, GSE121810, Gide2019, and Liu2019), (3) pre-treatment anti-PD-L1 (GSE165252, GSE176307, GSE183924, Mariathasan2018, and Miao2018), (4) on-treatment anti-PD-L1 (GSE165252), (5) pre-treatment anti-CTLA-4 (GSE115821, GSE140901, GSE165278, Gide2019, Miao2018, and VanAllen2015), and (6) on-treatment anti-CTLA-4 (GSE115821, GSE165278, and Gide2019) cohort. Due to the relatively small sample sizes, on-treatment datasets were excluded and three separate analyses were performed using the pre-treatment samples only. A summary of each dataset is listed in Table 1.

20Different technological, and clinical approaches were utilized in the involved studies such as RNA-sequencing (RNA-Seq) with Illumina or Ion Torrent platforms, and NanoString nCounter. Regarding therapies, many studies used combination therapy e.g., anti-PD-1, or anti-PD-L1, with or without anti-CTLA-4 therapy, or anti-PD-1 monotherapy, with or without anti-CTLA-4 treatment, while others investigated monotherapies only – such as ipilimumab, atezolizumab, or durvalumab.

21First, the pre-treatment samples in group #1 were investigated by computing ROC AUC and P-values for 29,755 genes. Following Bonferroni correction, values reaching more than P = 1.6E-06 were excluded from further analysis, which led us to 912 significant genes. The most significant hits included MARCKS (FC = 2.2, AUC = 0.724, P = 4.2E-12), SPIN1 (FC = 1.6, AUC = 0.682, P = 9.1E-12), PAK5 (FC = 1.8, AUC = 0.677, P = 6.1E-08), DDAH1 (FC = 1.7, AUC = 0.671, P = 3.4E-10), ELOVL6 (FC = 2.1, AUC = 0.667, P = 3.5E-10), SRC (FC = 1.6, AUC = 0.667, P = 5.9E-10), SETD7 (FC = 1.7, AUC = 0.663, P = 1.0E-09), FGFR3 (FC = 2.1, AUC = 0.657, P = 3.7E-09), YAP1 (FC = 1.6, AUC = 0.655, P = 6.0E-09), BACE1 (FC = 1.6, AUC = 0.654, P = 7.8E-09), STK35 (FC = 1.7, AUC = 0.651, P = 1.4E-08), TEAD3 (FC = 1.7, AUC = 0.649, P = 4.1E-08), TMPRSS4 (FC = 1.5, AUC = 0.640, P = 2.5E-07), PP1CB (FC = 1.6, AUC = 0.639, P = 2.0E-07), and BCL2 (FC = 2.2, AUC = 0.634, P = 9.7E-08) which were all upregulated in the non-responder group (Fig. 2). The complete list of all significant genes is provided in Supplementary Table S1. Notably, non-protein coding genes (such as pseudogenes, long intergenic non-protein coding RNAs, antisense RNAs, regulatory RNAs, small nucleolar RNAs, open reading frame, etc.) were excluded from our screening.

22By using all significant genes (n = 912) for GO analysis, mechanisms such as retrograde transport, vesicle recycling within Golgi (GO:0000301), ncRNA catabolic process (GO:0034661), and T cell receptor signaling pathway (GO:0050852) were significantly overrepresented among these genes (Supplementary Table S2).

23ROC AUC and P-values from 26,819 genes were computed and following Bonferroni-correction, values over P = 1.8E−06 were excluded. This way, we identified 38 significant genes. The complete list of all significant genes can be found in Supplementary Table S3. The Gene Ontology analysis shows that mechanisms connected to the C-type lectin receptor signaling pathway (GO:0002223), cellular response to lectin (GO:1990858), and positive regulation of natural killer cell-mediated cytotoxicity (GO:0045954) were overrepresented in the list of significant genes (n = 38) (Supplementary Table S2). There were no upregulated, druggable genes capable to predict resistance against anti-PD-L1 therapy.

24In this third analysis, ROC AUC and P-values were calculated for 22,561 genes in the pre-treatment group. Of these, 80 genes reached significance after Bonferroni correction. Among them, BLCAP (FC = 1.7, AUC = 0.735, P = 2.1E−06) was found as a druggable gene overexpressed amongst non-responding patients (Fig. 2). The complete gene list can be found in Supplementary Table S4. Non-protein coding genes were also excluded from this group. The GO analysis with multiple testing correction delivered no significant classification for these genes.

25We also investigated the power of established cancer biomarkers for predicting therapeutic response to immune checkpoint inhibitors. For this, a previously described panel of cancer biomarkers was utilized (see the complete list with results in Supplementary Table S5). The analysis was performed with two cohorts including pre-treatment anti-PD-1, and pre-treatment anti-CTLA-4 samples. Following Bonferroni-correction, only one gene, ALK (Anaplastic Lymphoma Receptor Tyrosine Kinase, FC = 1.6, AUC = 0.612, P = 3.3E−05) showed a correlation with anti-PD-1 response. In case of anti-CTLA-4, CD19 (FC = 1.6, AUC = 0.666, P = 6.5E−04) and PGR (Progesterone Receptor, FC = 1.2, AUC = 0.639, P = 6.8E−03) were found to be predictive. Notably, CD274 (PD-L1), as expected, was overexpressed in responding patients, however, failed to pass the significance threshold after Bonferroni-correction (P = 1.6E−06).

26Finally, we aimed to determine to what extent one can predict sensitivity to the anti-PD-1 dostarlimab using the integrated database of published datasets. For this, we performed ROC analysis to evaluate the anti-PD-1 biomarkers previously published in a cohort of rectal cancer. We used transcriptomic data of 419 samples from melanoma, bladder, renal cell, and gastric cancer in the anti-PD-1 pre-treatment cohort (n = 776). In this analysis, MLH1 and MSH6 achieved high predictive values (FC = 1.5, AUC = 0.682, P = 2.1E−11 and FC = 1.4, AUC = 0.629, P = 7.4E−06, respectively). Notably, 218 genes reached even higher ROC AUC than MLH1 (Fig. 3).

27Here, we integrated available data from multiple datasets and used this combined database to uncover biomarkers related to response against ICIs in three independent clinical settings. An advantage of the presented analysis pipeline is the utilization of real-world patient data. While most of the individual studies have only a limited number of samples, our combined patient cohort with well over a thousand patients provides a higher statistical power.

28Among the most significant genes related to resistance against anti-PD-1 treatment, we identified several potentially druggable targets. Fibroblast Growth Factor Receptor 3 (FGFR3) plays a pivotal role in tumorigenesis and the regulation of innate and adaptive immune systems. Overexpression of FGFR3 is associated with an immunologically cold, T cell-depleted phenotype, which contributes to a low ICI response rate in bladder cancer—just like a low PD-L1 expression in an FGFR3-mutant scenario. From multikinase inhibitors (such as anlotinib, dovitinib, lenvatinib, and ponatinib) to selective FGFR inhibitors (e.g., erdafitinib, infigratinib, pemigatinib), various small molecule inhibitors are available for solid tumors and lymphohematopoietic cancers. The combination of FGFR and PD-1 inhibition might also be beneficial. The Src Proto-Oncogene Non-Receptor Tyrosine Kinase (SRC) is a well-known oncogene contributing to cell growth, cell proliferation, and survival. Tumor-induced cytokines, e.g., Macrophage Inflammatory Protein 1, and 2 (MIP-1 and MIP-2) activate Src kinases in immune cells which lead to the production of pro-inflammatory cytokines (e.g., interleukin-1ß and 6, Tumor Necrosis Factor α) that activate cancer cells in a positive feedback loop. Multikinase inhibitors, such as the FDA-approved bosutinib, dasatinib, ponatinib, and vandetanib are currently being used for the treatment of chronic myeloid leukemia, acute lymphoblastic leukemia, and patients with thyroid cancer. B-Cell Lymphoma 2 (BCL2) is a major regulator of the “apoptotic-orchestra” by inhibiting apoptosis and promoting cell survival. Immune checkpoint molecules themselves promote an anti-apoptotic phenotype—leading to immune tolerance and low response rates. Venetoclax is the only FDA-approved small-molecule inhibitor against BCL2 in acute myeloid leukemia (AML), and chronic lymphatic leukemia—yet other drugs might follow both in hematologic and solid tumor malignancies. Yes-Associated Protein 1 (YAP1) is a transcriptional coactivator, which upon binding to many transcription factors, such as TEAD3 (Transcriptional Enhanced Associate Domain 3), regulates the Hippo-signalling pathway, contributing to tumor growth, and resistance. The Hippo-YAP system interferes with the innate immune system in multiple ways such as inhibiting the production of type I interferons (IFN-α, IFN-ß) and reactive oxygen species (ROS), attenuating NF-κB activation, or enhancing TNF-α and IL-1ß production. These events contribute to the suppression of the innate immune system thus escaping immune recognition – which eventually leads to tumor survival. Verteporfin (VP) is widely used for the treatment of macular degeneration, however, current studies highlighted the antitumor effects of VP either with photoactivation or without it. The capability of CA3, Super-TDU, statins, sitagliptin, SRC, FAK (Focal Adhesion Kinase), and tankyrase inhibitors to disrupt the YAP-TEAD complex were also discussed previously. SET Domain Containing Lysine Methyltransferase 7 (SETD7) has a broad target-specificity as it is involved in many biological processes by interacting with p53, Estrogen Receptor-Alpha (ERα), or YAP1. For this reason, SETD7 can both activate and inhibit tumor-survival signals. Upon methylation on K494 by SETD7, YAP1 accumulates in the cytoplasm and blocks the Hippo-pathway. This leads to the nuclear accumulation of ß-catenin and the activation of the Wnt/ß-catenin pathway which is one of the most well-known oncogenic pathways. Besides having a direct effect on cell proliferation and stemness, ß-catenin promotes a non-inflammatory milieu in tumors by inhibiting the activation and recruitment of CD8+ T cells and enhancing the infiltration and survival of regulatory T cells (Tregs). This leads to resistance to ICIs, emphasizing the potential of Wnt/ß-catenin as a predictive biomarker, or as a therapeutic target. Moreover, SETD7 can methylate p65 and thereby inhibit the expression of NF-κB, and is a positive regulator of Transforming Growth Factor Beta (TGF-ß) production – all these contribute to tumorigenesis. There are some promising results of inhibiting SETD7 with DC-S100, DC-S285, cyproheptadine, and ®-PFI-2. Spindlin 1 (SPIN1) is a histone methylation reader contributing to the epigenetic regulation of many oncogenic pathways so it is not surprising that SPIN1 was found to be overexpressed in many cancers. Notably, SPIN1 acts as a transcriptional coactivator of ß-catenin and T cell Factor 4 (TCF-4) enhancing their contribution to Wnt/TCF-4 pathways, which leads to tumor progression. Inhibitors of SPIN1 are being studied—e.g., A366, EML405, MS31, 4-aminoquinazoline and quinazolinethione derivatives, or VinSpinIn.

29We have found only one gene related to anti-CTLA-4 resistance: BLCAP or Bladder Cancer Associated Protein. As an apoptosis-inducing factor, BLCAP can initiate apoptosis in many tumors and is considered a tumor suppressor gene. Lost expression or degradation of BLCAP is observed in urothelial, renal, and cervical cancer, osteosarcoma, colorectal carcinoma, and human tongue carcinoma. Nonetheless, poor survival of bladder cancer patients correlates with strong nuclear expression of BLCAP, which is impacted by the interaction of BLCAP and Signal Transducer and Activator of Transcription 3 (STAT3) in the JAK/STAT-pathway. Selective pharmacological inhibition of BLCAP can be observed with aristolochic acid in vitro.

30In locally advanced rectal cancer, the loss of mismatch-repair genes was a highly significant biomarker of response to the anti-PD-1 inhibitor dostarlimab. When re-evaluating the previously published genes, we were able to confirm a significant correlation with survival in our cohort for MLH1 and MSH6. Nevertheless, 218 genes reached even higher significance than MLH1 in these patients, pointing out that other biomarkers might be even more suitable to select patient cohorts for immunotherapy. The samples in our patient cohort stemmed from melanoma, bladder, renal cell, and gastric cancer suggesting that the loss of mismatch-repair genes could also enhance sensitivity to anti-PD-1 therapy in these tumor types.

31There are some limitations of our study. Despite the importance of this topic, only a relatively small number of datasets have been found and included in the analysis. Since the currently used drugs are targeting proteins, an addition of protein-level data would have been useful. Unfortunately, we have not found even one dataset with additional protein abundance data. Finally, as ICIs are approved for advanced cancers, most patients have already received multiple treatment regimes, which might have resulted in significant background noise at the transcriptomic level preventing the identification of the most reliable biomarker candidates.

32In summary, we have established a database consisting of gene expression and clinical response data by combining 1434 solid tumor tissue samples obtained before or after immune-checkpoint inhibitor treatment. The most significantly upregulated, pharmacologically important (druggable) genes were identified in connection with the resistance against anti-PD-1, and anti-CTLA-4 treatments. Our extended analysis platform can help to identify, validate, and rank future biomarker candidates.