Checkpoint Genes at the Cancer Side of the Immunological Synapse in Bladder Cancer.

Immune checkpoint inhibitors have revolutionized cancer therapy, but not all cancers respond to the currently available drugs, and even within cancers considered responsive to such modality, response rates range between 15 and 40%, depending on the cancer type, the line of treatment, and yet unknown clinical/molecular factors. Coordinated expression of checkpoint proteins was shown to occur on T cells, probably allowing fine-tuning of the signal transmitted to the cell. We performed a bioinformatic analysis of the expression of putative checkpoint mRNAs at the cancer side of the immunological synapse from the bladder cancer tumorgenome atlas (TCGA) database. Fifteen mRNAs, corresponding to both coinhibitory and costimulatory checkpoints, were shown to be expressed above a designated threshold. Of these, seven mRNAs were found to be coexpressed: CD277, PD-1L, CD48, CD86, galectin-9, TNFRSF14 (HVEM), and CD40. The expression of 2 of these mRNAs-BTN3A1 (CD277) and TNFRSF14 (HVEM)-was positively correlated with overall survival in the TCGA database. All these seven mRNA share putative binding sites of a few transcription factors (TFs). Of these, the expression of the TF BACH-2 was positively correlated with the expression of checkpoint mRNAs from the network. This suggests a joint transcriptional regulation on the expression of checkpoint mRNAs at the bladder tumor side of the immunological synapse.

Introduction

There is an ongoing revolution in clinical oncology in the last decade following the realization that cancer develops an entire range of mechanisms to evade the host's immune response [1]. Extensive research is aimed at studying the cellular interface between cancer or antigen presenting cells (APCs) and lymphocytes, designated “the immunological synapse.” Immune checkpoint proteins—namely, transmembrane proteins coexpressed at both the cancer/APC and the lymphocyte side of immunological synapse—serve to modulate the signal transmitted from the cancer to the T cell, leading to either proliferation and activation (a costimulatory effect) or anergy and exhaustion (a coinhibitory effect) [2].

Three families of monoclonal antibodies targeting checkpoint inhibitors are already approved and being used to treat cancer—anti-CTLA4 (targeting the coinhibitory protein CTLA4 on T cells), anti-PD1 (targeting the coinhibitory protein PD-1 on T cells), and anti–PD-1L (targeting the coinhibitory protein PD-1L on cancer cells). Notwithstanding these major advancements, not all cancers respond to the currently available immune checkpoint inhibitors, and even within cancers considered responsive to such modality, response rates range between 15 and 40%, depending on the cancer type, the line of treatment, and yet unknown clinical/molecular factors.

Urothelial carcinoma of the bladder has long been perceived to be an immunogenic malignancy, and indeed intrabladder immune modulation with Bacillus Calmette–Guérin has been the mainstay of treatment for high-risk nonmuscle invasive bladder cancer (BLCA) for decades. Recently, both anti-PD1 and anti–PD-1L antibodies were shown to have activity in metastatic urothelial carcinoma of the bladder, with response rates ranging between 16 and 25%, depending on the trial and the agent [[3], [4], [5]]. There are now many ongoing clinical trials with combinations of immune checkpoint modulators in BLCA. Clearly, a better understanding of the immunological synapse in BLCA is warranted to advance immunotherapeutic treatment in this disease.

Previously, the expression of costimulatory and coinhibitory checkpoint proteins on the surface of T cells was shown to be coordinated and concerted (reviewed in Ref. [6]), symbolically metaphorized to resemble a tide wave of checkpoint protein activation [7]. Our hypothesis was that a similar coordinated expression of checkpoint proteins may occur on the cancer side of the immunological synapse, potentially allowing fine-tuning of the signal transmitted to the T cells. Here, we provide a bioinformatic analysis of coexpression networks of checkpoint mRNAs, based on the cancer genome atlas (TCGA) database,1 providing insights on potential new checkpoint genes that mandate further experimental research in this disease. We also point to a few transcription factors (TFs) that may be involved in the regulation of expression of these checkpoint genes, thus suggesting new potential targets for anticancer therapies.

Materials and Methods

Data Acquisition and Preprocessing

mRNA expression and metadata of the cases were acquired from The Cancer Genome Atlas (TCGA) [1] database using the “TCGAbiolinks” package in R [8]. We obtained the BLCA dataset (TCGA-BLCA), which consists of 412 tumor samples. For correlation aliases, we used the results of “HTSeq-FPKM” workflows for mRNA. We discarded duplicates of the same patient and selected only for experiments where matching sample was profiled for gene expression and miRNA expression. This resulted in a cohort of 405 tumor samples, which was later used for all correlation analyses. We also analyzed 21 samples of normal bladder urothelium (“controls”). Of note, healthy tissue samples marked as “controls” were derived from the same patients from which the tumor samples were derived.

Based on the literature, 22 genes potentially expressed at the cancer/APC side of the immunological synapse were selected for this analyses: PD-L1, PD-L2, CD40, CD48, CD70, CD80, CD86, CD112, CD137L, CD200, CD276, CD266, Gal-9, OX40L, HVEM, PVR, ICOSL, VTCN1, GITRL, VISTA, HHLA2, and IDO1. A comprehensive list of these checkpoints, their versatile names, and their suggested functions at the synapse is given in Table 1.

Table 1

A List of All Potential Checkpoints on the Cancer/APC Side if the Immunological Synapse.

On a Tumor	Additional Names	On an Immune System Cell (T Cells Mostly)	Coinhibitory (−) or Costimulatory (+) or Unknown (?)
MHC I/II		TCR	Signal 1
MHC I/II		KIR	−
MHC I/II		LAG3 (CD223)	−
CD80	B7-1, B7-H1, CD28LG, CD28LG1, LAB7, BB1	Binds to CD28 (+) or CTLA-4 (−)
CD86	B7-2, B7-H2, FUN-1, BU-63, B70, LAB72, CD28LG2,	Binds to CD28 (+) or CTLA-4 (−)
CD276	B7-H3, B7-3	Binds to CD28 (+) or CTLA-4 (−)
VTCN1	B7-H4, B7-4	Binds to CD28 (+) or CTLA-4 (−)
VISTA	B7-H5, B7-5, PD-1H, Gi24, VSIR	?
HHLA2	B7-H7	TMIGD2 (IGPR-1)	−
PD-L1	CD274, PDCD1LG1, B7-H1	Binds to PD-1	−
PD-L2	CD273, PDCD1LG2, B7-DC	Binds to PD-1	−
CD277	BTN3A, BT3.1, BTF5		?
CD134L	TNFSF4, OX40L	Binds to OX40 (TNFRSF4)	+
CD137L	TNFSF9, 4-1BBL	Binds to 4-1BB receptor (TNFRSF9)	+
CD70	TNFSF7, CD27L	TNFRSF7, CD27	+
B7RP1	ICOSL	ICOS	+
CD112	NECTIN2, HVEB, PVRL2	TIGIT	−
		CD226, DNAM-1	+
		CD112R	−
CD200	MOX1, MOX2, MRC, OX2	CD200R	−
CD48	BCM-1, BLAST-1	CD244, 2B4	+
Gal-9	LGALS9, galectin-9, LGALS9A, HUAT, Tumor Antigen HOM-HD-21	Binds to galactosides, HAVCR2, TIM3	−
GITRL	TNFSF18, TL6, HGITRL, AITRL, TNLG2A	TNFRSF18, GITR, AITR, CD357	+
CD40	TNFRSF5, P50	TNFSF5 (CD40LG, CD154, IGM, IMD3, HIGM1, T-BAM, GP39)	+
HVEM	TNFRSF14, HVEA, TR2, LIGHTR, ATAR, CD270	BTLA, HVEML, LIGHT, CD270L, TNFSF14	−
PVR	PVC, CD155, NECL5, TAGE4, HVED	VTN (vitronectin)
		CD96	−
		CD226, DNAM-1	+
		TIGIT	−
IDO	IDO1, TDO	↓ tryptophan	−

Correlation Analyses

We selected mRNAs of interest and created an n x m numeric matrix expression measures of reads per million per kilobase of transcript for RNA-seq or normalized qPCR values for PCR-based assays, where n is number of experiments (405 in case of BLCA cohort) and m is number of preselected mRNAs. We then generated a nonparametric correlation matrix using Spearman's rank-order correlation coefficient implemented in R. The statistical significance of each correlation is determined using a correlation test—a t-test is applied to the individual correlations using the following formula: t = r* sqrt(n-2)/sqrt(1-rˆ2). This method is implemented in “psych” package. The p values are then corrected to q-values using the false detection rate (FDR) method [9]. The resulting correlations and associated q-values are visualized for further inspection using “corrplot” package in R—our variation of the plot shows correlation values in lower triangle and q values in circles in the upper triangle of the heatmap-like plot [10].

Graphical Model Estimation

Correlation graphs represent the correlation matrix with nodes that indicate genes of interest and edges that represent correlation values. Green edges indicate positive correlations and red edges negative ones. The width of the edges and their color saturation corresponds to the absolute value of correlations and scale relative to the strongest weight in the graph. The graphs are organized as “spring” layout, which uses the Fruchterman–Reingold algorithm [11] to obtain a force-directed layout. In this solution, each node (connected and unconnected) repulse each other, and connected nodes also attract each other. After a number of iterations (500), a final logout is reached—the distance between the nodes corresponds well to correlation between the nodes—correlated nodes are close to each other, while anticorrelated (negative correlation) nodes are moved to distant parts of the graph.

Boxplots of Expression Values

Boxplots shows the FPKM values from primary solid tumor (PT) and normal tissue (NT). The p-values estimating statistical significance of difference between PT and NT are calculated using Mann–Whitney–Wilcoxon U-test.

Motif Analyses

Gene promoter sequences were obtained from Homo sapiens reference genome (hg38) using the Biostrings package in R. We defined promoters as regions 2000 bp upstream to 200 bp downstream annotated transcription start site (TSS) as customary (https://www.cell.com/molecular-cell/fulltext/S1097-2765(15)00262-2). Gene coordinates were obtained directly from TCGA annotations used in the RNA-seq pipeline. To find and analyze motifs, we used tools from MEME suit12: MEME—for de novo motif discovery [12], CentriMo [13]—to assess the central motif enrichment at promoter regions and TomTom—to find similar motifs in known motif databases. We ran two types of analyses—long motif search, where we de novo searched the whole 2.2 kb sequence at once, and short motif search, where we splatted the promoter sequence into 33 tiles of 100 bp that overlap by 25 bp. For each MEME output, we calculated central enrichment and annotated the motif using the search against known motif databases.

Results

The expression of 22 checkpoint mRNA suggested to be on the cancer/APC side of the immunological synapse is seen in Figure 1. Of these, 21 are transmembrane checkpoints (excluding IDO1 that codes for an intracellular protein). Fifteen had a median expression level of above 1 transcript per million, chosen here as the expression cutoff (Table 2). These include the coinhibitory molecules: PD-L1, CD276, galectin-9, CD200, and TNFRSF14 (HVEM); the costimulatory molecules: CD40, CD48, ICOSLG and CD137L; the molecules VTCN1, CD112, CD86, and PVR, probably implicated in both types of interactions (depending on their T-cell counterpart); and VISTA and CD277, the effect of which is currently uncertain. Of these, the median expression of CD276, PVR, and CD112 is significantly higher in tumor than in NT, and the median expression of VISTA, CD48, and CD200 is significantly lower in tumor than in NT (Figure 2).

Figure 1

Expression of checkpoint mRNAs in bladder cancer and normal samples in the TCGA database.

Table 2

Median Expression of 15 Checkpoint mRNAs in Bladder Tumor Samples in the TCGA Database.

mRNA Name	Median TPM
PVRL2	95.0
CD276	71.0
LGALS9	67.2
HVEM	49.5
CD40	20.4
VISTA	16.9
PVR	13.5
BTN3A1/CD277	12.9
CD48	6.3
ICOSLG	5.0
CD200	3.8
VTCN1	3.5
CD86	3.4
CD274/PD-1L	1.1
TNFSF9/CD137L	1.1

Figure 2

Differential analysis of 15 selected genes. Expression level of genes involved in the formation of the immunological synapse in bladder cancer tumor tissue and in normal bladder tissue in the TCGA database.

We then checked for correlations in the expression of these 15 mRNAs, calculating the Spearman rho correlation coefficient of each possible pair. Using a rho correlation coefficient cutoff of ≥0.5 to be statistically significant (corresponding to an adjusted q value of 0.1 or less following FDR correction for multiple comparisons), 7 mRNAs were found to be coexpressed: CD277, PD-1L, CD48, CD86, galectin-9, TNFRSF14 (HVEM), and CD40 (Figure 3). The coexpression network is depicted in Figure 4. To verify that these seven genes are indeed significantly coexpressed, we checked all 6435 combinations of 7 genes out of 15. The connectivity of each of these 6435 potential networks was scored using 3 methods—calculating the mean correlation of the network; using the network coefficients described by Onnela and colleagues; and using the network coefficients of Zhang and Horvath (the last two described in Ref. [14]) According to all 3 methods, the abovementioned seven genes were significantly more connected that any other combination, with mean p values of 0.00062, 0.0047, and 0.018 for the 3 methods, respectively.

Figure 3

Correlation matrix with Spearman rho values (bottom triangle) and q-values (upper triangle, in circles) of 7 checkpoint mRNAs, based on TCGA data. Tonation and size of the values are associated with the strength of correlation.

Figure 4

Expression correlation network of selected 7 checkpoint mRNAs with the rho values ≥ 0.5, based on Spearman rho coefficients. Green line indicates positive correlations, and line thickness reflects the strength of the correlation.

Of these 7 checkpoint mRNAs, the expression of 2 of these mRNAs—BTN3A1 (CD277) and TNFRSF14 (HVEM)—was correlated with overall survival in the TCGA database (Figure 5).

Figure 5

Kaplan–Meier survival curves of bladder cancer patients from the TCGA with mRNA expression of HVEM (left) or CD277 (right panel) above (red) or below (green) the median.

Last, we searched the putative promoter sites for potential TF binding motifs common to all of the seven coexpressed genes. Because in most cases the promoter is not known, the standard promoter length used for analyses was 2000 base pairs (bp) upstream and 200 bp downstream of the TSS. The following TFs were found to have putative motifs in the promotors of all seven coexpressed mRNAs: MAFK, BACH-1, BACH2, NFE2L2, TFAP2A, TFAP2B, TFAP2C, MEIS2, MEIS3, PKNOX1, PKNOX2, TGIF1, and TGIF2 (Figure 6). To assess whether these 13 motifs are indeed significantly overrepresented in these 7 promotors, we performed a hypergeometric test for the overrepresentation of these motifs against a set of 1000 randomly selected promoter sets of 7 genes. We found that all 13 motifs were significantly overrepresented, with p-values ranging from 3.51*10−5 for MEIS1 to 0.026 for TFAP2B. This indicates that the occurrence of these common motifs in the 7 putative promotors cannot be explained by mere change and suggests that indeed the TFs take part in regulating the coexpression of these 7 checkpoint genes.

Figure 6

Common TF-binding motifs in the putative promoters of the 7 coexpressed checkpoint mRNA with their orientation and calculated p value.

Of these 13 TFs, only the expression of BACH-2 was positively correlated with the expression of 2 checkpoint mRNAs with the network —namely CD86 and CD48—using the same rho cutoff of 0.5. The expression of BACH2 was also correlated with the expression of the CD200, which was in itself borderline correlated with the expression of CD48. There was a nonsignificant trend of positive correlation between the expression of BACH-2 and PD-L1 (CD274; rho = 0.4, q = 0.2). Figure 7 depicts a comprehensive network of all checkpoint mRNAs and BACH-2 mRNA.

Figure 7

Expression correlation network of all checkpoint mRNAs (red nodes) and the BACH-2 transcription factor (blue node), based on Spearman rho coefficients. Green lines represent positive correlations, and red lines represent negative ones. Line thickness represents the strength of the correlation.

Discussion and Conclusions

Using bioinformatic tools and analyses, we show here that of 21 transmembrane checkpoint genes potentially expressed on the cancer side of the synapse, the mRNA of 15 is expressed above a designated threshold in BLCA samples from the TCGA. These include genes with potentially costimulatory and coinhibitory functions. Of these, the expression of 7 is highly correlated to one another, and a high expression of 2 of these—HVEM and CD277—is associated with better prognosis in the TCGA database. All the putative promotor sites of these 7 mRNAs contain potential binding motifs for several TFs. Of these TFs, the expression of BACH-2 has a positive correlation with the expression of several components of the checkpoint mRNA network.

Our analysis reveals that the mRNAs CD276, galectin-9, HVEM, PVRL2, and PVR, all suggested to have coinhibitory or mixed functions, are expressed at significantly higher levels than PD-1L mRNA. There have been dozens of clinical trials with agents targeting the PD1/PD-1L axis in recent years, and several of these monoclonal agents have been approved for treatment of BLCA in the last year or two. Still, most patients do not respond to single agent anti-PD1/PD-1L inhibition in BLCA, mandating a deeper understanding of the intricate immunological synapse in this disease.

A recent review summarized the expanding repertoire of potential targets for immune modulation in the bladder and pointed out CD276, HHLA2, and VTCN1 as such [15]. Our results corroborate this review as well as previous research showing the CD276 mRNA and protein are expressed in BLCA [16] and contribute to the accumulating evidence of its role in BLCA. CD276 expression promotes invasion and metastasis of urothelial cancer cells in vitro [17], and clinical trials with the anti-CD276 antibody MGA271 [18] are under way (www.clinicaltrials.gov). HHLA2 mRNA was not expressed in the TCGA database of BLCA, and VTCN1 was only expressed at low levels. Galectin-9 protein expression was also shown to be associated with prognosis in BLCA [19]. HVEM, PVRL2, and PVR have not been implicated in BLCA so far to our knowledge, but recently encouraging clinical activity of a monoclonal antibody against nectin-4 (PVRL4) conjugated to a microtubule-disrupting agent was presented [20]. Our results further strengthen the notion that PD-1L is not necessarily the most important or potent coinhibitory checkpoint in BLCA.

The expression of seven of these mRNAs—CD277, CD48, CD86, PD-1L, galectin-9, HVEM, and CD40—was found to be significantly correlated to one another. This coexpression resembles the coexpression of checkpoint mRNAs at the T-cell interface of the synapse [6,21]. If these results are corroborated at the protein level, then their implication is that coinhibition occurs concurrently by several checkpoints. In such circumstances, targeting a single checkpoint coinhibitor, or even two, may not suffice, as is indeed seen in most patients.

Three TFs have putative binding sites in the promotors of the coexpressed checkpoint genes: MAFK, BACH2, and NFE2L2. BACH2 TF has recently been reported as crucial in several pathways related to the immune system functioning [22,23]. It is also involved in NF-κB signaling pathway that is in itself cardinal in BLCA pathogenesis [23]. Aberrant function of BACH2 has been already implicated in several cancer types, such as lymphomas [24]. MAFK is a powerful transcription regulator, in a form of leucine zipper, cooperating with NFE2 TFs family, hence acting as a transcription silencer or enhancer, depending on the protein interactions [24,25]. Interestingly, MAFK may act as a competitive repressor of some NFE2 TFs, blocking their binding sites [[25], [26], [27]].

The correlation between BACH-2 expression and the checkpoint mRNA network, taken with its putative ability to regulate the expression of these mRNAs, warrants further experimental research. If indeed checkpoint mRNAs, and consequently proteins, are coregulated, then this may open novel avenues for modulating the expression of checkpoints at the immunological synapse.

The expression of two mRNAs—HVEM (TNFRSF14) and CD277 (BTN3A1)—was positively associated with survival of BLCA patients. CD277 is a member of the butyrophilin subfamily 3 that shares significant sequence similarities and predicted common structural features with other members of the B7 superfamily. It was shown to be expressed in the microenvironment of ovarian cancer [28] and pancreatic adenocarcinoma [29]. Butyrophilin 3 A isoforms are critical activating molecules of Vg9Vd2 T cells [29], and their association with improved outcome in urothelial cancer may indirectly hint to the importance of Vg9Vd2 T cells in eliciting an antineoplastic response in this tumor type.

HVEM (TNFRSF14) interacts with multiple ligands expressed in the immune system including the TNF superfamily cytokines, and its interactions with its ligands can be altered in pathologic settings, resulting in dysregulated immune response. It has mainly been implicated in several types of lymphoma and also in nonhematological cancers such as melanoma. Its function as a coinhibitor or costimulator of T cells seems to be context-dependent (reviewed in Ref. [30]), and a specific role for it in BLCA has not yet been described.

There are several important limitations to our work. First, our results are based on retrospective analysis of a single, albeit large and well-curated, database. Our proposed networks should therefore be verified by others as well. Second, our analysis is based on the expression of mRNAs and not proteins. A comprehensive analysis of protein expression and protein–protein interactions at the cancer cell membrane is mandatory. Second, the mRNA expression measured in the TCGA samples does not necessarily stem from expression on the cancer cells but is rather a composite measurement of expression on all cell types within the microenvironment; having said that, we claim that it still represents a summation of tumor checkpoint mRNA expression as a whole. Last, we used here arbitrary, if customary, cutoffs for both expression and coexpression. Using different cutoffs would have slightly altered the components of the suggested networks, but not the essence of our findings.

With the advancement of the field of cancer immunotherapy, the main target of treatment has shifted from cancer cells themselves to entire cancer microenvironments, including immune cells and stromal components. Addressing the particular tumor characteristics, changing cancer immune set points and modulating the delicate interactions between tumor cells and immune cells are currently all major strategies aimed at enhancing the patient's ability to combat the disease (reviewed, for example, in Ref. [31]). Our work contributes to the body of accumulating evidence as to the nature of tumor–immune–cell interactions at the immunological synapse and opens up novel avenues for further experimental and translational research.

Funding

This work was partly funded by a grant from the Israeli Cancer Association (ICA) and from a generous donation from the Jack Craps foundation.