Mapping Immune Correlates and Surfaceome Genes in BRAF Mutated Colorectal Cancers

Simple Summary The identification of novel therapeutic strategies for Colorectal Cancer (CRC) patients with BRAF mutations is mandatory, since most of the current treatments against this tumor type, including novel compounds, show limited efficacy. In this article, we describe upregulated proteins in the surface of cells within the tumor that can be used as targets to specifically guide novel treatments. In addition, we also observed that the transcriptomic profile matches with antigen presenting cells, such as dendritic cells and macrophages having “antigen processing and presentation of exogenous peptide antigen via MHC class II” as main molecular function, favoring an immunoreactive microenvironment. Therefore, the combination of anti PD(L)1, together with other co-inhibitor receptors from the ones presented in this manuscript, should be explored to treat BRAF mutated CRC patients. Abstract Despite the impressive results obtained with immunotherapy in several cancer types, a significant fraction of patients remains unresponsive to these treatments. In colorectal cancer (CRC), B-RafV600 mutations have been identified in 8–15% of the patients. In this work we interrogated a public dataset to explore the surfaceome of these tumors and found that several genes, such as GP2, CLDN18, AQP5, TM4SF4, NTSR1, VNN1, and CD109, were upregulated. By performing gene set enrichment analysis, we also identified a striking upregulation of genes (CD74, LAG3, HLA-DQB1, HLA-DRB5, HLA-DMA, HLA-DMB, HLA-DPB1, HLA-DRA, HLA-DOA, FCGR2B, HLA-DQA1, HLA-DRB1, and HLA-DPA1) associated with antigen processing and presentation via MHC class II. Likewise, we found a strong correlation between PD1 and PD(L)1 expression and the presence of genes encoding for proteins involved in antigen presentation such as CD74, HLA-DPA1, and LAG3. Furthermore, a similar association was observed for the presence of dendritic cells and macrophages. Finally, a low but positive relationship was observed between tumor mutational burden and neoantigen load. Our findings support the idea that a therapeutic strategy based on the targeting of PD(L)1 together with other receptors also involved in immuno-modulation, such as LAG3, could help to improve current treatments against BRAF-mutated CRC tumors.


Introduction
The identification within many tumors of druggable oncogenic vulnerabilities settled the basis for the development of a new generation of selective therapeutic agents that are currently in use in the clinical setting [1,2]. One of the strategies that has been used to act on these targets is based on the use of small molecule inhibitors that can permeate cell membranes and inhibit the function of intracellular proteins. Conversely, monoclonal neutralizing antibodies act by blocking soluble factors or plasma membrane proteins containing extracellular domains [3,4]. Examples of successful targeted therapies in cancer include, among many others, the small molecule inhibitors of the BRAF V600 oncogene in melanoma, and monoclonal antibodies against the Epidermal Growth Factor Receptor (EGFR) in non-small cell lung cancer (NSCLC) [1]. Another more recent successful example is the use of inhibitors that can selectively target the K-RAS p.G12C mutation in NSCLC [5,6]. The progressive implementation of these types of selective therapeutic strategies has revolutionized cancer treatment, leading to the development of more individualized therapies that have demonstrated significant improvements in the outcome of patients that harbor specific mutations [4].
Although Colorectal Cancer (CRC) is one of the most prevalent tumors, the identification of druggable molecular alterations lags behind other less frequent solid tumors, such as melanoma or NSCLC, where the use of selective therapies has clearly demonstrated clinical efficacy [7]. Recently, several druggable mutations have been reported in genes frequently mutated in CRC, including those affecting BRAF or K-RAS, among others [8].
In the case of K-RAS mutations, their presence in CRC predicted a lack of response to anti-EGFR antibodies [9].
In CRC, B-RafV600 mutations have been identified in 8-15% of the patients and were associated with detrimental outcomes and a lack of response to anti-EGFR inhibitors [10]. Treatment with BRAF inhibitors such as Encorafenib has shown activity, particularly if combined with anti-EGFR antibodies [11]. However, although this and other similar combinational therapeutic strategies are promising, only a minor fraction of patients exhibited positive responses to these treatments. Moreover, most of the patients who respond to this therapeutic approach have a limited extension in progression-free survival [10,11]. In addition to this limitation, the toxicity profile of these agents could in some cases impair their administration [12]. In this context, the identification of novel druggable opportunities or therapeutic options to improve activity and reduce toxicity is mandatory.
Lately, immunotherapy has gained momentum by demonstrating impressive clinical responses in various cancer types [13]. Several strategies have been assayed to influence or reactivate the immune system of cancer patients in a manner that can contribute to fight malignant recurrent diseases. One of these strategies is interfering with the inhibitory signal that blocks the host immune response on cancer cells by using anti PD(L)1 antibodies. Another option is the use of antibodies against membrane receptors to boost antibody drug dependent cytotoxicity (ADCC) or complement dependent cytotoxicity (CDC) [14]. An additional approach is the use of antibodies against membrane proteins to vectorize small molecule inhibitors or chemotherapeutic agents as antibody drug conjugates (ADCs) [15,16]. With this approach, off-target non-tumor toxicity will be considerably reduced as the therapy is directly guided to the cancer cell [17].
The surfaceome is considered the set of all plasma membrane proteins that have extracellular domains. In a cancer context, the clarification of the surfaceome of the cell types that form the tumor mass, including cancer and stromal cells, is of paramount interest as it can help to identify specific therapeutic targets and therefore facilitate the design of novel individualized therapies [18].
In this work, we took advantage of transcriptomic data to explore BRAF-mutated CRC. Following this strategy, we identified a surfaceome signature that could help to predict the immune status of these tumors and therefore contribute to design more effective therapies against this CRC subtype.

Identification of BRAF Mutations in CRC Patients, Data Collection and Processing
We used data contained at cBioportal (www.cbioportal.org) (accessed on 1 February 2022) [19,20] (TCGA dataset) to explore all mutations in the BRAF gene in patients with colorectal cancer. This web resource also provides mutated variants mapped to genomic domains. Protein expression in the cell membrane was identified using the Human Surfaceome Atlas (https://wlab.ethz.ch/surfaceome/) (accessed on 15 March 2022) [21].

Functional Annotation of De-Regulated Genes
We used the publicly available EnrichR online platform (https://maayanlab.cloud/ Enrichr/) (accessed on 16 March 2022) [22] to address the gene ontology biological process and molecular function related to each gene set. We represented the most relevant pathways according to their adjusted p-value.

Outcome Analysis
The KM Plotter Online tool [23] (https://kmplot.com/analysis/, last accessed on 20 March 2022), was used to evaluate the relationship between up-regulated genes' expression and clinical outcome in patients with rectum adenocarcinoma. This open access database contains 165 samples and allowed us to investigate Free Progression (FP) and Overall Survival (OS) of up-regulated genes in the rectum adenocarcinoma subtype. False discovery rate (FDR) indicates replicable associations across multiple studies.

Expression Analysis
The analysis comparing the expression level of individual genes in BRAF mutated samples compared with wildtype ones was carried out using data from the Cancer Dependency Map (DepMap) portal (https://depmap.org/portal/, accessed on 28 October 2021) for cell lines and the Xena Functional Genomics Explorer web server from the University of California Santa Cruz (https://xenabrowser.net/, last accessed on 25 March 2022) [24] for TCGA colon and rectum adenocarcinoma samples.

Correlation between Gene Expression and Immune Cell Infiltration
To explore the associations between gene expression and immune infiltration cells we used TIMER2.0 (http://timer.cistrome.org/, accessed on 5 April 2022). TIMER provides 4 modules (Gene, Mutation, sCNA and Outcome) to explore the association between immune infiltrates and genomic changes [25]. The gene-correlation module was used to link gene expression with activation of T cell markers.

Correlation with Tumor Mutational Burden and Antigen Load
The total number of Mutect2-identified somatic mutations was used to compute the "mutational burden" in all TCGA colon and rectal cancer patients [26]. Neoantigen load was determined using the NetMHCpan-3.0 [27], as described previously [28]. Finally, Spearman rank correlation coefficients were computed to correlate mutational burden and neoantigen load and gene expression.

Graphical Design
Bars, dot plots, heatmap and volcano plots were represented using GraphPad Prism software (GraphPad Software 9.0.0, San Diego, CA, USA).

Mapping Transcriptomic Differences in CRC Tumors Based on BRAF Mutated Status
To identify upregulated genes characteristically expressed in CRC BRAF mutated tumors, we interrogated open datasets as described in the materials and methods section. Of 396 patients, 62 of them harbored BRAF mutations, corresponding to 15.65% of all patients with CRC. Using a threshold fold change (FC) equal or higher than two, we selected 210 genes that were upregulated and 280 genes that were downregulated when comparing BRAF mutated CRC patients versus BRAF wt ones. Among them, 57 genes were upregulated and presented at the surface of the cells according to the Surfaceome Atlas [21].
We displayed in a volcano plot those genes up and downregulated in BRAF mutated versus BRAF wt CRC tumors ( Figure 1A). We observed that downregulated genes were more abundant, being 57.14% of total dysregulated genes, compared to 42.86% of upregulated genes. However, the FC in upregulated genes was higher, having a mean fold change of 3.35 versus 3.16 in the downregulated ones. When focusing only on those highly upregulated, we identified GP2  Upregulated genes that were expressed at the surfaceome are presented in Figure 1C. We noticed that they had a mean FC of 2.23. This group included the previously described genes, in addition to NTSR1 (5.66 FC), VNN1 (4.94 FC), and CD109 (4.52 FC) ( Figure 1D). Supplementary Figure S1 displays the expression of the mentioned genes in CRC patients with wt and mutated BRAF. In CRC, cell lines AQP5, CLDN18, CD109, and NTSR1 were highly expressed in the BRAF mutated group (Supplementary Figure S2).

Functional Analysis of Identified Upregulated Genes
Next, we evaluated the function of the listed upregulated genes. As can be seen in Figure 2A, the main biological functions included the following: "antigen processing and presentation of exogenous peptide antigen as biological function" and "antigen processing and presentation of exogenous peptide antigen via MHC class II". Genes within these three functions included CD74, LAG3, HLA-DQB1, HLA-DRB5, HLA-DMA, HLA-DMB, HLA-DPB1, HLA-DRA, HLA-DOA, FCGR2B, HLA-DQA1, HLA-DRB1, and HLA-DPA1 (Figure 2A). When evaluating the molecular functions, the most prevalent functions included MHC class II receptor activity and MHC class II protein complex binding, including similar genes ( Figure 2B). We noticed that all these genes express proteins that have a strong protein-protein interaction (PPI) (with a PPI enrichment p-value < 1.0 × 10 −16 ), forming a cluster of 13 nodes with 64 edges (when the expected number was one), meaning they have an average node degree of 9.85 with an average local clustering coefficient of 0.938 ( Figure 2C). This means that these proteins have more interactions among themselves than what would be expected for a random set of proteins. Therefore, such enrichment indicates that the proteins are at least partially biologically connected, as a group.   Figures S3 and S4). The statistically positive correlation (Rho > 0.5, p-value < 0.005) between the selected surfaceome genes and other immune populations is displayed in Supplementary Figure S5, where a positive correlation with macrophages was also observed in the colon adenocarcinoma (COAD) dataset of the rectum. Of note, antigen presentation is also produced by macrophages. At the same time, this set of genes presented a negative correlation with myeloid derived suppressor cells (MSDCs, see Supplementary Figure S5). Finally, we observed that the mentioned CD74 signature correlated with MSI-H tumors and the colon cancer subtype CMS1 that represents 14% of the tumors and is associated with immune strong activation ( Supplementary Figures S6 and S7). Next, we aimed to correlate the presence of the identified genes with the expression of PD(L)1 or its receptor PD1. As shown in Figure 4, we observed a strong correlation within the presence of the selected genes, particularly for CD74, HLA-DPA1, and LAG3 with PD1 and PD(L)1 in both COAD and READ datasets. These data suggest that an immunosuppressive microenvironment is present with a high expression of exhausted T cells (PD1 positive). Of note, the genes that were highly upregulated and described in Figure 1 were not associated with PD (L)1 (Supplementary Figure S8) or any specific immune cell subtype (Supplementary Figure S9). These findings, in addition to the purity score, suggest that these genes correspond to the tumor cells and are not involved in the modulation of the immune system.

Correlation with Tumor Mutational Burden and Antigen Load
Then, we integrated the genes by using their mean expression into two signatures termed CD74 and AQP5 signatures and correlated these signatures with the mutation burden (Supplementary Figure S10). In this analysis, both signatures had a low but positive correlation (correlation coefficients (or corr. coeff.) = 0.193, p = 3.9 × 10 −5 and corr. coeff. = 0.23, p = 1.1 × 10 −6 , for the CD74 and AQP5 signatures, respectively). The same analysis was also performed for the neoantigen load by using the number of expressed peptides binding to MHC molecules, and both signatures reached a similar level of significance (corr. coeff. = 0.267, p = 4.1 × 10 −8 , and corr. coeff. = 0.244, p = 5.3 × 10 −7 for the CD74 and AQP5 signatures, respectively).

Discussion
In the present article, we explored surfaceome genes in CRC patients that harbor mutations in the BRAF gene with the main goal to identify potential immunologic druggable opportunities.
BRAF is an oncogene found to be mutated in a significant proportion of tumors: up to 10-15% of all CRC patients [8]. The BRAF V600 is the most frequent mutation, representing 98% of the mutations in the BRAF gene [8]. This mutation is mainly observed in the right sided CRC, and in poorly differentiated tumors [29]. In addition, tumors with this mutation primarily belong to the consensus molecular subgroup (CMS) CMS1 or microsatellite instability (MSI) immune subtype that integrates microsatellite instability high (MSI H) tumors and associates with an immunologic transcriptomic profile [30]. However, although 42% of BRAF mutated CRC tumors can be included within this group, there is still a high proportion that belong to other subtypes not characterized as immune-related [31].
In our study, we observe that CRC tumors with BRAF mutations express a high number of genes linked with antigen presenting functions. All of them belonged to the major histocompatibility complex (MHC) II molecules that are mainly expressed in professional antigen-presenting cells. Antigens presented in these molecules are extra-cellular proteins that are endocytosed and digested in the lysosomes, and those peptide fragments are exposed by the MHC II molecules [32,33]. Of note, high antigen load correlates with an adaptive immune response and, therefore, with a higher efficacy of immune checkpoint inhibitors to be active [34,35]. Conversely, a high neoantigen load does not guarantee a good response to immune checkpoint inhibitors (ICIs) due to the presence of an immune suppressive microenvironment [35,36].
In line with the previous finding, when matching the immune populations with the transcriptomic profile of BRAF-mutated CRC patients we observed that the best correlation was with antigen presenting cells, such as dendritic cells and macrophages. This finding indicates that BRAF mutated cells favor a microenvironment enriched in antigen presenting cells. In addition, we also observed a moderate correlation with CD8 and a minor correlation with CD4 T cells, and a strong association with PD1 and PD(L)1 expression. The mentioned CD74 signature, that is associated with colon tumors with MSI-H and with the colon cancer subtype termed CMS1 (that represents 14% of all tumors), has been characterized as MSI Immune and it is known to be hypermutated and microsatellite unstable [30]. These findings suggest that, although there is a presence of effector T cells in the tumors, these cells are inhibited or exhausted. Furthermore, some of the genes identified as upregulated in BRAF-mutated CRC patients were associated with detrimental outcomes, including HLA-DOA, CD74, HLA-DMB, and HLA-DQA1, further supporting the idea that the microenvironment of these tumors contributes to the development of an aggressive immunosuppressive phenotype. We therefore postulate the idea that novel strategies aimed at improving the efficacy of ICI in BRAF-mutated CRC patients should combine the use of anti-PD(L)1 antibodies with inhibitors of other receptors involved in the regulation of immune responses that could help to revert this immunosuppressive phenotype. In this context, it is worth underlining the relevance of LAG3 as a negative regulator of effector T cells as well as of antigen presenting cells [37,38]. Interestingly, an anti-LAG3 antibody termed Relatlimab was recently demonstrated to be active in combination with nivolumab as a first line treatment of melanoma [39]. It is therefore tempting to speculate that a similar strategy could be undertaken in BRAF mutated CRC.
In addition to the changes observed in the expression of the immune system-associated genes described above, we also found a strong upregulation of other genes encoding for proteins that contain domains that are exposed at the outer layer of the plasma membraneincluding GP2, CLDN18, AQP5, TM4SF4, NTSR1, VNN1, and CD109. The expression of these genes does not correlate with immune populations or PD(L)1 expression, suggesting that their expression is restricted to tumor cells. This finding could also be of great interest since it opens the possibility of using some of these proteins as therapeutic targets against BRAF-mutated CRC. For example, antibodies against CLDN18 (encoding for the tight junction protein Claudin 18) are currently under development in gastric and pancreatic cancer [40]. Our observations now set the basis for the potential evaluation of these antibodies alone or in combination with anti-PD(L)-1 or anti-LAG-3 antibodies in CRC tumors harboring BRAF mutations. Likewise, VNN1 has been linked with detrimental overall survival [41] and could also be also therapeutically explored.
On the other hand, we acknowledge that our study has limitations. The first one relates to the lack of discrimination between the stromal and tumor compartment. In this context, single cell sequencing could differentiate between the populations where these genes are expressed. Finally, we consider that data presented here should be validated in the laboratory or even in a more accurate manner by using immunohistochemistry (IHC) analysis in human samples.

Conclusions
Taken together, our findings support the notion that BRAF-mutated colorectal tumors favor a strong immune activated state enriched with antigen-presenting cells and T lymphocytes, and therefore associates with MSI-H and CSM1 CRC tumors. We hypothesize that a therapeutic strategy based on the blockade of PD(L)1, as well as other receptors also involved in immuno-modulation, such as LAG3, together with the pharmacological inhibition of BRAF, could help to improve current treatments against BRAF-mutated CRC tumors.

Institutional Review Board Statement:
Given the fact that our study did not involve the use of human samples, as it was a bioinformatic study, no ethical approval was required.
Informed Consent Statement: Patient consent was waived due to the fact that no human samples were used in the study: Public available bioinformatic data was used.
Data Availability Statement: All the data used in this manuscript is open and would be shared upon request to the corresponding authors.