Novel approach to identify putative Epstein-Barr-virus microRNAs regulating host cell genes with relevance in tumor biology and immunology.

The human Epstein-Barr virus is associated with several human solid and hematopoietic malignancies. However, the underlying molecular mechanisms including virus-encoded microRNAs (miRs), which lead to the malignant transformation of infected cells and immune evasion of EBV-associated tumors, have not yet been characterized. The expression levels of numerous known EBV-specific miRs and their suitability as diagnostic and/or prognostic markers were determined in different human EBV-positive tissues followed by in silico analyses to identify putative EBV-miR-regulated target genes, thereby offering a suitable screening strategy to overcome the limited available data sets of EBV-miRs and their targeted gene networks. Analysis of microarray data sets from healthy human B cells and malignant-transformed EBV-positive B cells of patients with Burkitt's lymphoma revealed statistically significant (p < 0.05) deregulated genes with known functions in oncogenic properties, immune escape and anti-tumoral immune responses. Alignments of in vivo and in silico data resulted in the prediction of putative candidate EBV-miRs and their target genes. Thus, a combinatorial approach of bioinformatics, transcriptomics and in situ expression analyses is a promising tool for the identification of EBV-miRs and their potential targets as well as their eligibility as markers for EBV detection in different EBV-associated human tissue.

Introduction

The human Epstein–Barr virus (EBV) is a double-stranded DNA virus belonging to the γ-subfamily of the human herpes viruses (HHV), which is known to exert a certain oncogenic potential.[1] Indeed, EBV infection is associated with several human malignancies including B cell Non-Hodgkin lymphoma (NHL), like Burkitt’s lymphoma (BL)[2] or NK/T-cell lymphoma,[3] classical Hodgkin’s lymphoma (cHL),[4] gastric carcinoma[5] and nasopharynx carcinoma (NPC).[6]

The viral gene expression pattern differs between latent and lytic phase of the EBV infection, which could be further subgrouped based on their gene expression profile.[7] In detail, in the latency 0 phase present in infected and circulating memory B cells in vivo, only non-protein coding viral transcripts are produced. Latency 1 phase typically observed in BL exerts EBNA1 expression along with the non-coding viral transcripts, while the latency 2 phase characteristic for cHL and NPC exhibits an additional expression of LMP1 and LMP2. In contrast, in the latency 3 phase of B cell transformation, EBNA2, EBNA-LP, BHRF1, the three latent membrane proteins LMP1, 2A and 2B in addition to the non-coding transcripts can be detected.[7]

The viral genome encodes for >80 genes and to date for 48 viral microRNAs (miRs), which are classified according to their genomic location (BamHI rightward reading frame 1 (BHRF1) and BamHI-A rightward transcript (BART)) into ebv-miR-BHRFs and ebv-miR-BARTs.[8-10] During the distinct EBV infection phases, a differential expression pattern of EBV-miRs was described.[7,11] MiRs are small RNAs with approximately 22 nucleotides (nts) in length and regulate sequence-specific the post-transcriptional gene expression by binding preferentially, but not exclusively to the 3′ untranslated region (UTR) of mRNAs leading to translational inhibition and/or degradation of the targeted mRNA.[12,13] Concerning canonical miR target/mRNA interactions featuring complementarity to the seed sequence, the second to eighth nt of the miR is critical for miR targeting, but recent studies also reported miR-target mRNA interactions without perfect seed pairing thereby possibly limiting the reliability of the in silico prediction tools.[14] Most of the EBV-miRs are involved in the development of EBV-driven malignancies and exert a variety of important regulatory functions in EBV-mediated tumorigenesis and cancer progression.[8,15] Screening for EBV-miR targets is even more complex, since the miR target prediction tools available do not include the sequences of the EBV-miRs.

Due to the regulatory potential of tumor-relevant target genes, the expression of certain EBV-miRs has been shown to be of prognostic relevance in EBV-associated malignancies. For example, in gastric cancers, the ebv-miR-BART20-5p is linked to a poor patients’ survival,[16] while the ebv-miR-BART7 has been described as potential plasma marker in NPC patients.[17]

The implementation of miRs as additional prognostic markers in formalin-fixed and paraffin-embedded (FFPE) tumor tissues of EBV-associated neoplasia is of growing importance due to their increased stability when compared to mRNAs or transcripts >80 nucleotides (nts) in length, which rapidly degrade in FFPE tissue specimens. Indeed, the quality and stability of miRs extracted from FFPE specimen was comparable to that of miRs derived from fresh tissues and was not affected by different fixation times.[18-20]

Depending on the functions of the targeted mRNAs, miRs can exert tumor suppressive or even oncogenic functions.[21] Furthermore, the miR-dependent regulation of several (tumor) immunologic relevant target genes has even led to the definition of immune modulatory miRs involved in immune recognition and/or immune escape of tumors.[22] Furthermore, EBV-encoded molecules, like miRs (in particular miR-BART7-3p) and proteins (LMP1), are directly involved in the epithelial to mesenchymal transition (EMT), which in solid tumors might allow virus-infected cells to migrate into other tissues and locations in combination with a deregulated immune surveillance.[23-26]

Recently, virus-encoded miRs have been shown to be involved in the various interactions between the EBV and the human host cells, thereby affecting different cell biological processes including immune surveillance, apoptosis induction, cell proliferation, tumor growth/suppression and cell migration.[27-32]

The exact identification of the molecular targets and signaling pathways accompanying these processes is of general importance for the understanding of the EBV-driven malignant transformation and immune escape. Therefore, this study aims to identify deregulated relevant target genes with tumor biological importance, which reflect direct or indirect targets of EBV-encoded miRs thereby highlighting the combinatorial approach of in silico analyses with transcriptome analyses and in situ expression data.

Materials and methods

Cell culture, RNA extraction, cDNA synthesis and qPCR

The cell lines were obtained from the ATCC (www.atcc.org) and cultured in RPMI1640 medium (Gibco™, Fisher Scientific, Waltham, MA, USA) supplemented with 10% (V/V) fetal bovine serum (Anprotec, Germany), 1 x penicillin, streptomycin (Sigma-Aldrich, Taufkirchen, Germany), 1 x sodium pyruvate (Gibco™, Fisher Scientific, Germany), 1 x L-glutamine (Lonza, Basel, Switzerland) and 25 mM (final concentration) Hepes buffer (c-c-pro, Oberdorla, Germany) at 37°C and 5% CO2. RNAs of the EBV-positive cell lines were kindly provided by the Prof. Ofer Mandelboim laboratory.

Total RNA was extracted with TRIzol™ (Invitrogen™, Fisher Scientific, Waltham, MA, USA) and reverse transcribed by RevertAid First Strand cDNA Synthesis Kit (Fisher Scientific) according to Chen and co-authors[33] and Jasinski-Bergner and co-authors[12] followed by qPCR analyses as already described in detail.[12] For the determination of the relative miR expression values, the human endogenous hsa-miR-3960 was applied as housekeeping gene. All primers used and their conditions are listed in Table 1.

Table 1.

List of the oligonucleotides used for the identification of EBV-miRs

Primer	Application	Sequence (5′→3′)	Conditions
stemloop ebv-miR-BART3-5p	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACAGCACA	42°C
stemloop ebv-miR-BART4-3p	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACACACCT	42°C
stemloop ebv-miR-BART7-3p	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACCCCTGG	42°C
stemloop ebv-miR-BART8-3p	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACTCTACG	42°C
stemloop ebv-miR-BART10-3p	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACACAGCC	42°C
stemloop ebv-miR-BART13-3p	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACTCAGCC	42°C
stemloop ebv-miR-BART15	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACTCAAGG	42°C
stemloop ebv-miR-BART17-5p	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACCTTGTA	42°C
stemloop ebv-miR- BHRF1-3	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACTGTGCT	42°C
stemloop hsa-miR-3960	cDNA synthesis	GTCGTATCCAGTGCAGGGTCCGAGGTATTCGCACTGGATACGACCCCCCG	42°C
qPCRfw ebv-miR-BART3-5p	qPCR	GCCCACCTAGTGTTAGTGT	60°C
qPCRfw ebv-miR-BART4-3p	qPCR	GCCCCGCACCACTAGTCACC	60°C
qPCRfw ebv-miR-BART7-3p	qPCR	GCCCCATCATAGTCCAGTGT	60°C
qPCRfw ebv-miR-BART8-3p	qPCR	GCCCGTCACAATCTATGGGGT	60°C
qPCRfw ebv-miR-BART10-3p	qPCR	GCCCTACATAACCATGGAGTT	60°C
qPCRfw ebv-miR-BART13-3p	qPCR	GCCCTGTAACTTGCCAGGGAC	60°C
qPCRfw ebv-miR-BART15	qPCR	GCCCGTCAGTGGTTTTGTTT	60°C
qPCRfw ebv-miR-BART17-5p	qPCR	GCCCTAAGAGGACGCAGGCA	60°C
qPCRfw ebv-miR- BHRF1-3	qPCR	GCCCTAACGGGAAGTGTGTA	60°C
qPCRfw hsa-miR-3960	qPCR	GCCCGGCGGCGGCGGAGGC	60°C
qPCR pan rev	qPCR	GTGCAGGGTCCGAGGT	60°C

fw: forward; rev: reverse.

Patient samples, EBER-CISH staining, immunohistochemistry and miR extraction from FFPE specimen

FFPE specimen from EBV-positive B cell lymphoma (n = 7), T cell lymphoma (n = 2), cHL (n = 2), NPCs (n = 12) and tonsils of patients with infectious mononucleosis (n = 13) as well as the respective EBV-negative controls (EBV-negative lymph nodes, n = 9; EBV-negative NPC adjacent mucosa tissues, n = 4; EBV-negative tonsils, n = 7) were provided as anonymised slides from the archive of the Institute of Pathology, Martin-Luther-University Halle-Wittenberg, Halle, Germany and conducted according to the principles expressed in the Declaration of Helsinki.

Scientific research on FFPE tissue samples from patients was approved by the Ethical Committee of the Medical Faculty of the Martin Luther University Halle-Wittenberg, Halle, Germany (2017–81).

The Epstein–Barr virus-encoded RNA chromogenic in situ hybridization (EBER-CISH) was performed on 4 µm thick FFPE tissue slides using the ZytoFast PLUS CISH Implementation Kit AP-Permanent Red (ZytoVision, Bremerhaven, Germany) according to the manufacturer’s instructions. Immunohistochemistry (IHC) was performed on a Bond III automated immunostainer (Leica Biosystems Nussloch GmbH, Wetzlar, Germany) using the Bond Polymer Refine Detection Kit (DS9800-CN). The primary antibody directed against LMP1 (clone CS.1–4, Agilent, Santa Clara, CA, USA) was applied as recommended by the manufacturer.

From each FFPE tissue block, three to four consecutive tissue slides (each 4 μm thick) were pooled for further miR extraction. Therefore, the MasterPure™ Complete DNA and RNA Purification Kit (Lucigen, Middleton, USA) was utilized according to the manufacturer’s protocol. The cDNA synthesis and qPCR analyses were performed as mentioned above.

In silico analyses of microarray data

The free online R2: Genomics Analysis and Visualization Platform (http://r2.amc.nl) was employed to determine the gene expression of various marker genes. Two data sets, both based on the same microarray chips (u133p2, Affymetrix Genechips, Santa Clara, CA, USA), were used: (i) healthy human B cells sorted for CD19[34] and (ii) the EBV-positive B cell line P493-6 derived from a human BL patient[35] cultured under hypoxic conditions, which contains a tetracycline-inducible myc transgene.[36]

After selection of the data set[34,35] at the R2 database, the option “view a gene” was selected. An Analysis was performed by submitting the respective gene symbols to obtain the exact expression levels in every single study sample.

Statistical evaluation and software

All statistically significant (two-tailed Mann-Whitney-U test) values were marked with * (p < 0.05), with ** (p < 0.005) and with *** (p < 0.0005). For statistical analyses and visualization Microsoft Excel 2016 (Microsoft Corporation, Redmond, WA, USA) and Prism GraphPad9 (GraphPad Software, San Diego, CA, USA) were applied.

Results

Identification of EBV-miR expression patterns in human EBV-positive tissue specimen of distinct origins and their role as diagnostic biomarkers

To identify host cell genes, which are targeted by EBV-encoded miRs expressed during latency, a selection of nine different EBV-miRs known to be expressed during latency was analyzed in three EBV-positive B cell lines (DAUDI, EB1, RAJI) in comparison to three EBV-negative B cell lines (BV173, CA46, DG75) by qPCR. The relative copy numbers were determined using the very abundantly expressed endogenous hsa-miR-3960 as housekeeping gene. While the EBV-encoded miRs were not detectable in the EBV-negative B cell lines (data not shown), all EBV-miRs were expressed in the EBV-positive cell lines. Some of the EBV-miRs exerted a heterogeneous expression within the three EBV-positive cell lines, while others were identified in all three EBV-positive cell lines with the highest expression of ebv-miR-BART7-3p and ebv-miR-BART10-3p (Figure 1A). To validate the postulated suitability of EBV-miRs as putative diagnostic markers for EBV detection in FFPE specimen as alternative and/or in addition to the EBER-CISH staining, the expression of the selected nine EBV-miRs was quantified by qPCR in a panel of EBV-positive human FFPE specimen. The EBV status of these FFPE specimen was proven by EBER-CISH staining and LMP1 expression determined by IHC as exemplarily shown in Figure 1B and 1C, respectively, and summarized in Table 2.

Figure 1.

Expression and quantification of selected EBV-encoded factors in vitro and in vivo. (A) qPCR-based quantification of selected EBV-miRs in three EBV-positive cell lines. The relative copy numbers of the selected EBV-miRs were determined as described in Materials and Methods. The human endogenous hsa-miR-3960 served as housekeeping gene for calculation of the relative expression values, which are expressed as mean of three biological replicates. (B) representative EBER-CISH staining of an EBV-positive tonsil derived from a patient with infectious mononucleosis (sample ID: IM3). EBER-CISH staining was performed according to Materials and Methods. (C) representative anti-LMP1 immunohistochemistry staining of an EBV-positive tonsil derived from a patient with infectious mononucleosis (sample ID: IM3). LMP1 expression was determined by IHC using an anti-LMP1 antibody as described in Material and Methods.

Table 2.

Comparison of different EBV detection methods using EBV-positive malignant tissues

Sample ID	Tissues analyzed	LMP1 positivity	EBER-CISH positivity	Detection of EBV-miRs
L1	lymphoma	yes	yes	yes
L2	lymphoma	no	yes	yes
L3	lymphoma	yes	yes	yes
L4	lymphoma	yes	yes	yes
L5	lymphoma	no	yes	yes
L6	lymphoma	yes	yes	yes
L7	lymphoma	yes	yes	yes
L8	lymphoma	yes	yes	yes
L9	lymphoma	no	yes	yes
L10	lymphoma	yes	yes	yes
IM1	infectious mononucleosis tonsil	yes	yes	yes
IM2	infectious mononucleosis tonsil	yes	yes	yes
IM3	infectious mononucleosis tonsil	yes	yes	yes
IM4	infectious mononucleosis tonsil	yes	yes	yes
IM5	infectious mononucleosis tonsil	yes	yes	yes
IM6	infectious mononucleosis tonsil	yes	yes	yes
IM7	infectious mononucleosis tonsil	yes	yes	yes
NPC1	nasopharyngeal carcinoma	yes	yes	yes
NPC2	nasopharyngeal carcinoma	no	yes	yes
NPC3	nasopharyngeal carcinoma	yes	yes	yes
NPC4	nasopharyngeal carcinoma	yes	yes	yes
NPC5	nasopharyngeal carcinoma	yes	yes	yes
NPC6	nasopharyngeal carcinoma	yes	yes	yes
NPC7	nasopharyngeal carcinoma	yes	yes	yes
NPC8	nasopharyngeal carcinoma	yes	yes	yes

The conventional detection methods for EBV infection using LMP1-specific immunohistochemistry and EBER-CISH were employed in different EBV-associated biopsies and compared to the qPCR-mediated detection of EBV-miRs.

In addition, EBV-miRs were analyzed by qPCR in the same specimen from EBV-positive B cell lymphoma, T cell lymphoma, cHL, NPC and tonsils of patients with infectious mononucleosis as well as respective EBV-negative controls (for further information see Materials and Methods section). As shown in Table 2, EBV-miRs were detected in all EBV-positive samples.

In detail, the expression of ebv-miR-BART7-3p, ebv-miR-BART10-3p, ebv-miR-BART17-5p and ebv-miR-BART15 was statistically significant higher when compared to the respective EBV-negative controls (Figure 2A–I). With the exception of ebv-miR-BHRF 1–3, the other investigated EBV-miRs showed increased expression levels in at least one of the analyzed tissues. However, the expression levels of the different EBV-miRs varied within the analyzed EBV-positive specimen of distinct origin suggesting a tissue-specific expression. Interestingly, a statistically significant difference was found for ebv-miR-BART4-3p (p = 0.0344) and ebv-miR-BART17-5p (p = 0.0026) between EBV-positive NPC and EBV-positive lymphoma samples. Despite this heterogeneity, these results support the hypothesis of the eligibility of certain EBV-miRs as markers for EBV detection in addition to EBER-CISH.

Figure 2.

Quantification of selected EBV-miRs in vivo in different EBV-positive and EBV-negative tissue samples. (A-I) qPCR-based quantification of selected EBV-miRs in EBV-positive tonsils of patients with IM, EBV-positive NPC tumors and EBV-positive lymph nodes of B cell lymphoma patients as well as respective EBV-negative controls was performed by determination of the relative copy numbers using the human endogenous hsa-miR-3960 as housekeeping gene as described in Materials and Methods. For calculation of statistically significance two-tailed Mann-Whitney-U test was performed.

Development of a strategy for identification of EBV-mediated deregulated, tumor-relevant genes

To investigate pathways, which might be altered during the process of EBV-induced malignant transformation, in silico analyses were performed. Since the standard in silico miR target prediction databases do not comprise EBV-miRs, including the candidate EBV-miRs, an alternative in silico screening strategy was established to identify potential human tumor relevant target genes, which are altered by EBV-encoded miRs expressed during latency.

Two microarray data sets available at the R2 microarray database (http://r2.amc.nl) were selected for comparison. The data set of Jima et al.[34] utilized HG-U133A gene chips (Affymetrix) for the analysis of flow cytometry sorted CD19-positive B cells of healthy human donors, while the study of Kim and co-authors[35] applied the same gene chips and investigated EBV-positive BL.[37,38]

Implementing in silico analysis of these model systems, the molecular impact of EBV-driven malignant transformation affecting different tumor-relevant pathways was determined to get insights into the regulation of malignancy-related genes in the generation of EBV-associated B cell neoplasia, but perhaps also in general for EBV-linked diseases including NPCs.

Analysis of components of the innate immune system and immune surveillance

First, the expression of genes relevant for the functionality of the innate immune system known to exert anti-viral activities were compared. These include involving anti-viral cytokines/chemokines, their receptors, downstream signaling molecules and pattern recognition receptors with anti-viral functions like the toll-like receptors (TLR)3, TLR7, TLR8, and TLR9 and the intracellular pattern recognition factors RIG-1 (DDX58), MDA5 (IFIH1) and LGP2 (DHX58)[39-42] (Figure 3A–B). A statistically significant enhanced expression was detected for the anti-viral chemokine CXCL10 (p = 2.1E-05), an almost statistically significant gene expression was detected for CXCL9 (p = 0.0703) and for the receptor CXCR3 (p = 0.0447) in EBV-positive BL cells (Figure 3C). In contrast, a statistically significant downregulation in BL cells was found for some components of the type I interferon (IFN) signaling, such as IFNAR1 (p = 0.0002), JAK1 (p = 0.0004) and STAT2 (p = 0.0011), while RIG-1 (DDX58; p = 0.0001) and MDA5 (IFIH1; p = 0.0018) were downregulated in BL cells (Figure 3C).

Figure 3.

Comparative analysis of the gene expression of anti-viral factors of the innate immune system derived from microarray data. Microarray data were analyzed as described in Material and Methods. The mean expression values were calculated (A) in healthy B cells and (B) in EBV-positive BL cells. Statistically significant deregulated genes were marked by a red arrow and almost statistically significant deregulated genes were marked by a blue arrow. These results were summarized in (C) including the calculated p values.

In addition, the expression of genes involved in the immune surveillance as well as immune evasion of tumors was analyzed (Figure 4A–C). These genes included the components of the antigen processing and presentation machinery (APM) and the NKG2D ligands MIC-A/B and UL16-binding proteins (ULBPs), which are frequently downregulated in tumors. It is noteworthy that the applied gene chips lack probes for ULBP5 and ULBP6. Furthermore, tumor-induced immune response-relevant molecules were determined. Interestingly, only two known immune modulatory genes, namely PD-L1 (CD274) and PVR (CD155 or NECL5), were statistically significantly upregulated (p = 0.0155 and p = 0.0384) in BL when compared to healthy B cells (Figure 4C).

Figure 4.

Comparative analysis of the gene expression patterns of molecules involved in immune surveillance and immune modulation derived from microarray data. Microarray data were analyzed as described in Material and Methods. The mean expression values were calculated for (A-C) genes involved in immune surveillance and in immune modulation. Statistically significant deregulated genes were marked by a red arrow. These results were summarized in (C) including the calculated p values.

Analysis of marker genes with prognostic relevance

EBV-positive BL cells exerted a statistically significant induced expression of some proliferation marker genes, like MKI67, E2F1 and BUB1 with the strongest effect for the cell cycle regulators MCM2 (p = 0.0091), MCM3 (p = 0.0033) and MCM4 (p = 0.0052) (Figure 5C).

Figure 5.

Comparative analysis of the gene expression patterns of proliferation markers of known oncogenes derived from microarray data. Microarray data were analyzed as described in Material and Methods. The mean expression values were calculated for (A-C) proliferation markers and for (D-F) oncogenes. Statistically significant deregulated genes were marked by a red arrow and almost statistically significant dysregulated genes were marked by a blue arrow. These results were summarized in (C) including the calculated p values.

Since oncogenes are frequently involved in the regulation of cell proliferation and cell growth, the gene expression of 14 common oncogenes was analyzed in the EBV-positive BL cells (Figure 5D–F). As shown in Figure 5D–F, a statistically significant upregulation of ABL (p = 2.1E-07), ERBB2 (p = 5.5E-06), NRAS (p = 2.7E-06), SYK (p = 2.4E 03) and an almost statistically significant upregulation of PDGFB (0.0579) were found in BL cells when compared to healthy B cells (Figure 5D–F).

Based on the assumption that metastasis formation is a major cause of death from cancer,[43] markers of the epithelial–mesenchymal transition (EMT) and the mesenchymal–epithelial transition (MET) were analyzed. This is based on the fact that EBV infection is not only linked to BL or cHL,[37,44] but also to solid tumor entities including inter alia NPCs and gastric cancers (GCs).[45-49] Since no data sets for solid EBV-associated tumor entities were available, the EMT and MET marker genes were also analyzed within the B cell and BL data sets, despite these results might only exert a surrogate for EBV-associated solid tumors.

It is known that at diagnosis some of the EBV-associated lymphoma are already at stage IV and infiltrated into the bone marrow, liver and lungs with differences within the lymphoma subtypes.[50] However, healthy B cells and the EBV-positive BL cells revealed a statistically significant downregulation of the two epithelial marker genes KLF4 (p = 0.0274) and KRT18 (p = 0.0338) (Figure 6A–B), while a statistically significant expression was found for FGFR1 (p = 4.7E-05), MMP14 (p = 0.0041), TGFB1 (p = 0.015), VIM (p = 0.0001) and ZEB2 (p = 0.0017) in healthy B cells (Figure 7A–C) and a borderline statistically significant higher gene expression for CD44 (p = 0.059), FOXC2 (p = 0.0657) and MMP2 (p = 0.0817) in BL cells. This borderline statistical significance at the transcript level might be considered for putative protein analyses.

Figure 6.

Comparative analysis of the gene expression from epithelial markers derived from microarray data. Microarray data were analyzed as described in Material and Methods. The mean expression values were calculated for (A-B) EMT markers. Statistically significant deregulated genes were marked by a red arrow. These results were summarized in (C) including the calculated p values.

Figure 7.

Comparative analysis of the gene expression from mesenchymal markers derived from microarray data. Microarray data were analyzed as described in Material and Methods. The mean expression values were calculated for (A-B) MET markers. Statistically significant dysregulated genes were marked by a red arrow and almost statistically significant deregulated genes were marked by a blue arrow. These results were summarized in (C) including the calculated p values.

Identification of EBV-regulated genes involved in the tumor metabolism

An altered metabolism in malignancies is a well-studied hallmark of cancer.[51] The increased glycolysis with a high catalysis and secretion of lactate in malignant tissues has been described as Warburg effect[52] contributing to an acidic tumor microenvironment (TME), thereby leading to a decreased activity of tumor-infiltrating immune effector cells and by acting as radical quencher enhancing the survival of tumor cells against radiation therapy and certain chemotherapies.[53] During EBV immortalization of B cells into LCLs, increased energetic demands for continued proliferation drive an upregulation of glycolysis and extracellular acidification.[54,55] Therefore, the gene expression of enzymes involved in the glycolysis pathway leading to the synthesis of lactate were compared between healthy B cells and BL cells (Figure 8A).

Figure 8.

Comparative analysis of the gene expression pattern for glycolysis-relevant enzymes focusing on to the Warburg effect from microarray data. Microarray data were analyzed as described in Material and Methods. The mean expression values were calculated (A) in healthy B cells and in EBV-positive BL cells. Statistically significant deregulated genes were marked by a red arrow and almost statistically significant dysregulated genes were marked by a blue arrow. These results were summarized in (B-C) including the calculated p values.

Indeed, 13 out of 19 investigated genes of glycolysis-relevant enzymes exerted a statistically significant increased expression and two genes showed an almost statistically significant gene expression in BL cells (Figure 8B–C). The strongest induction of glycolysis-relevant enzymes in BL cells was found for pyruvate kinase (p = 6.2E-06), triosephosphate isomerase (p = 0.0007), phosphofructokinase (p = 0.0001), glucose-6-phosphate isomerase (p = 0.0007) and hexokinase (p = 0.0001), which is accompanied by an increased proliferation of BL cells.

Chromosomal mapping of the deregulated tumor relevant target genes

In total, the expression of 161 cancer-related genes was analyzed between healthy B cells and malignant transformed EBV-positive BL cells. Compared to healthy B cells, statistically significant or almost significant altered gene expression patterns with either upregulated (n = 48) or downregulated (n = 5) genes by the EBV-driven malignant transformation of BL cells were found and summarized in Table 3 including their genomic localization. The deregulated genes were randomly scattered over the human chromosomes, but most genes were located on chromosome 12 (n = 6), chromosome 9 (n = 5) and chromosome 10 (n = 5). These results might be of importance due to known genetic alterations on these chromosomes within tumor diseases in general, which will be addressed in more detail within the discussion section. None of the deregulated genes were located on the chromosomes 5, 7, 18 and Y.

Table 3.

Identification of differentially expressed genes upon EBV-driven malignant transformation

Upregulated genes after malignant transformation						Downregulated genes after malignant transformation
Gene symbol		Function		Genomic localization		Gene symbol		Function		Genomic localization
CXCL10*		anti-viral chemokine		4q21.1		IFNAR1		receptor for anti-viral IFN-α		21q22.11
CXCL9		anti-viral chemokine		4q21.1		JAK1		component of the IFN downstream signaling		1p31.3
CXCR3		receptor for CXCL10 and CXCL9		Xq13.1		STAT2		component of the IFN downstream signaling		12q13.3
DDX58* (RIG-1)		anti-viral intracellular pattern recognition receptor		9p21.1		KLF4		epithelial marker		9q31.2
IFIH1* (MDA5)		anti-viral intracellular pattern recognition receptor		2q24.2		KRT18		epithelial marker		12q13.13
CD274* (PD-L1)		ligand for the inhibitory receptor PD-1 on T cells		9p24.1		Upregulated genes after malignant transformation
PVR* (CD155)		ligand for the inhibitory receptor TIGIT on NK cells		19q13.31		Gene symbol		Function		Genomic localization
MKI67		proliferation marker		10q26.2		HK1*		glycolysis relevant gene		10q22.1
E2F1*		proliferation marker		20q11.22		HK2*		glycolysis relevant gene		2p12
MCM2*		proliferation marker		3q21.3		GPI*		glycolysis relevant gene		19q13.11
MCM3*		proliferation marker		6p12.2		PFKL*		glycolysis relevant gene		21q22.3
MCM4*		proliferation marker		8q11.21		PFKM		glycolysis relevant gene		12q13.11
BUB1		proliferation marker		2q13		PFKP		glycolysis relevant gene		10p15.2
CD44*		mesenchymal marker		11p13		ALDOA*		glycolysis relevant gene		16p11.2
FGFR1*		mesenchymal marker		8p11.23		TPI1		glycolysis relevant gene		12p13.31
FOXC2*		mesenchymal marker		16q24.1		PGK1		glycolysis relevant gene		Xq21.1
MMP2*		mesenchymal marker		16q12.2		PGAM1		glycolysis relevant gene		10q24.1
MMP14*		mesenchymal marker		14q11.2		PGAM5		glycolysis relevant gene		12q24.33
TGFB1*		mesenchymal marker		19q13.2		ENO1*		glycolysis relevant gene		1p36.23
VIM*		mesenchymal marker		10p13		PKM		glycolysis relevant gene		15q23
ZEB2		mesenchymal marker		2q22.3		LDHA		glycolysis relevant gene		11p15.1
NRAS		oncogene		1p13.2		LDHB*		glycolysis relevant gene		12p12.1
PDGFB		oncogene/ growth factor		22q13.1		ABL1*		oncogene		9q34.12
SYK		oncogene		9q22.2		ERBB2		oncogene		17q12

*genes, which are reported in literature to be dysregulated upon EBV infection. Statistically significant or almost statistically significant upregulated (A) or downregulated (B) genes associated with EBV transformation are listed.

Implementation of a combinatorial strategy for prediction of novel putative EBV-miR-regulated target genes

Due to the lack of in silico EBV-miR prediction databases, the RNA sequences including the 5’-untranslated region (UTR), the coding sequence (CDS) and the 3’-UTR of 5 down-regulated cancer-related host genes were analyzed for putative EBV-miR binding sites, involving all EBV-miR sequences, currently available within the mirbase online database (www.mirbase.org)[56,57] by using the RNAhybrid prediction software (bibiserv.cebitec.uni-bielefeld.de/rnahybrid).[58,59] This software allowed the calculation of the free energy of the EBV-miR mRNA target interaction (Table 4), but these in silico calculations have some limitations, since they do not consider any specifics of viral miRs, possibly containing putative viral specific RNA modifications with relevance for mRNA binding affinity and targeting.[60] Five EBV-miRs exerted strong interaction sites to even more than two of the five down-regulated tumor-relevant mRNA targets (IFNAR1, JAK1, STAT2, KLF4, KRT18). ebv-miR-BART4-5p and −13-5p could per se be involved in the down regulation of all five target genes, while the EBV-miRs ebv-miR-BART6-3p and ebv-miR-BART17-3p could be regulatory for four and ebv-miR-BART7-3p for three target genes (Table 4). The ebv-miR-BART7-3p was even among the most abundantly expressed EBV-miRs analyzed in the EBV-positive BL cell lines as well as in the EBV-positive FFPE specimen (Figure 1A and Figure 2A) suggesting a strong contribution of this EBV-miR in malignant transformation.

Table 4.

In silico predicted free energies of putative EBV-miR-mRNA target interactions using RNA hybrid.[58,59]

miRNA	Sequence	Free energy of the in silico predicted miR-mRNA interaction [kcal/mol/]
		IFNAR1	JAK1	STAT2	KLF4	KRT18
EBV-miR-BART1-5p	ucuuaguggaagugacgugcugug	−30,5	−27,5	−30,9	−26,5	−29,8
EBV-miR-BART1-3p	uagcaccgcuauccacuauguc	−24,6	−28	−25	−27,6	−30,3
EBV-miR-BART2-5p	uauuuucugcauucgcccuugc	−24,8	−26	−24,5	−28,4	−23,8
EBV-miR-BART2-3p	aaggagcgauuuggagaaaauaaa	−22,5	−27,7	−22,8	−27,5	−25,3
EBV-miR-BART3-5p	accuaguguuaguguugugcu	−24	−25,3	−26,1	−24,4	−27,8
EBV-miR-BART3-3p	cgcaccacuagucaccaggugu	−28,6	−27,1	−28,1	−29,7	−27,3
EBV-miR-BART4-5p	gaccugaugcugcuggugugcu	−30,3	−31,8	−33,8	−31,5	−33,9
EBV-miR-BART4-3p	cacaucacguaggcaccaggugu	−30,3	−25,7	−28	−26,6	−25,1
EBV-miR-BART5-5p	caaggugaauauagcugcccaucg	−27,8	−25,8	−30,5	−29	−30,4
EBV-miR-BART5-3p	gugggccgcuguucaccu	−29,1	−30	−27,2	−32	−29,9
EBV-miR-BART6-5p	uaagguugguccaauccauagg	−24,3	−24	−25,9	−24,1	−23,8
EBV-miR-BART6-3p	cggggaucggacuagccuuaga	−31	−29,4	−30,6	−34	−31,8
EBV-miR-BART7-5p	ccuggaccuugacuaugaaaca	−26,3	−23,8	−23,8	−24,9	−26,9
EBV-miR-BART7-3p	caucauaguccaguguccaggg	−29,6	−30	−34,7	−31	−28,7
EBV-miR-BART8-5p	uacgguuuccuagauuguacag	−24,4	−24,2	−21,3	−23,8	−25,1
EBV-miR-BART8-3p	gucacaaucuauggggucguaga	−27	−28,7	−27,6	−28,6	−25,7
EBV-miR-BART9-5p	uacuggacccugaauuggaaac	−29,4	−25,4	−23,1	−25,1	−25,4
EBV-miR-BART9-3p	uaacacuucaugggucccguagu	−27,5	−27,9	−27	−25,9	−26,2
EBV-miR-BART10-5p	gccaccucuuugguucuguaca	−27,2	−28,2	−28,9	−29,3	−26,3
EBV-miR-BART10-3p	uacauaaccauggaguuggcugu	−27,8	−25,3	−25,4	−27,9	−26,4
EBV-miR-BART11-5p	ucagacaguuuggugcgcuaguug	−29,4	−27,6	−27,8	−27,4	−29,8
EBV-miR-BART11-3p	acgcacaccaggcugacugcc	−32,6	−28,6	−26,6	−31,4	−26,3
EBV-miR-BART12	uccugugguguuuggugugguu	−28,2	−29,9	−27,3	−31	−28,5
EBV-miR-BART13-5p	aaccggcucguggcucguacag	−31,1	−32,5	−30,7	−34,7	−30,5
EBV-miR-BART13-3p	uguaacuugccagggacggcuga	−29,2	−29,5	−35,1	−34,5	−29,9
EBV-miR-BART14-5p	uacccuacgcugccgauuuaca	−21,9	−26	−26	−27,8	−26,4
EBV-miR-BART14-3p	uaaaugcugcaguaguagggau	−25,6	−27	−32	−27,5	−24,6
EBV-miR-BART15	gucagugguuuuguuuccuuga	−26	−27,3	−25,3	−26,9	−24,6
EBV-miR-BART16	uuagauagagugggugugugcucu	−29,4	−30,2	−27,2	−27,8	−28,5
EBV-miR-BART17-5p	uaagaggacgcaggcauacaag	−26,5	−24,1	−23,9	−25,8	−24,6
EBV-miR-BART17-3p	uguaugccugguguccccuuagu	−30,8	−32,1	−34,3	−29,7	−30,2
EBV-miR-BART18-5p	ucaaguucgcacuuccuauaca	−22,9	−21,5	−23	−25,6	−20,3
EBV-miR-BART18-3p	uaucggaaguuugggcuucguc	−26,9	−29	−26,1	−27,8	−25,6
EBV-miR-BART19-5p	acauuccccgcaaacaugacaug	−23,9	−21,2	−22	−27,3	−21,9
EBV-miR-BART19-3p	uuuuguuugcuugggaaugcu	−24,5	−25,7	−23,2	−25,5	−25,4
EBV-miR-BART20-5p	uagcaggcaugucuucauucc	−25,4	−25,2	−32,3	−24,3	−24,6
EBV-miR-BART20-3p	caugaaggcacagccuguuacc	−26,3	−25,8	−25,6	−32,6	−23,9
EBV-miR-BART21-5p	ucacuagugaaggcaacuaac	−20,5	−21,1	−22,5	−18,2	−19,7
EBV-miR-BART21-3p	cuaguugugcccacugguguuu	−26	−28,1	−27,4	−27,1	−26
EBV-miR-BART22	uuacaaagucauggucuaguagu	−24,2	−26,8	−25,3	−22,3	−22,4
EBV-miR-BHRF1-1	uaaccugaucagccccggaguu	−30,6	−28,3	−29,3	−31,6	−28,4
EBV-miR-BHRF1-2-5p	aaauucuguugcagcagauagc	−25,8	−22,5	−24,8	−23,9	−22,2
EBV-miR-BHRF1-2-3p	uaucuuuugcggcagaaauuga	−28,1	−21,7	−22,8	−21,6	−20,6
EBV-miR-BHRF1-3	uaacgggaaguguguaagcaca	−26,8	−25,3	−25,4	−24,5	−22,6

The predicted free binding energy of 5 selected EBV-miR targets identified is given. In bold, free energies >30 are marked.

In summary, this screening strategy allows to identify the most abundantly expressed EBV-miRs in the context of deregulated tumor biologic-/tumor immunologic-relevant target genes, which have to be subjected to additional molecular biological and functional studies to prove their predicted correlation and even the direct interactions.

Discussion

A putative host cell mRNA directly bound by viral encoded miRs should lead to a decrease of the respective protein and/or to the decay of the miR-bound mRNA. In addition, upregulated host cell gene expression could be caused by indirect effects of EBV-miRs, e.g. by direct binding and downregulation of suppressors/negative regulators of the respective gene. Based on this assumption, a downregulation of gene expression can be caused by miRs indirectly, e.g. by direct binding and downregulation of an inducer/positive regulator. In order to identify genes regulated by EBV-miRs, the expression of EBV-miRs was monitored in EBV-positive and EBV-negative tissues of distinct origin demonstrating even tissue-specific differences.

The highest expressed EBV-miRs were identified and were used as a prerequisite to reduce the number of putative EBV-miR candidates for follow up in silico analyses.

Next to the EBV-miRs, also the EBV-encoded lnc RNAs EBERs and BARTs are known to interfere with the gene expression of infected host cells. Recent studies report a shift towards an oncogenic phenotype by EBER1 and EBER2 overexpression in in vitro cell lines, but also in in vivo murine models,[61-63] which was due to an enhanced AKT signaling and an increased c-MYC (MYC) expression. In our study, an upregulated MYC expression was found in BL cells when compared to healthy B cells, which was not statistically significant (p = 0.57). The EBV-encoded protein LMP1 could also activate the PI3K-AKT signaling pathway thereby interfering with DNA repair mechanisms.[64] These observations strengthened the involvement of additional molecules in the process of virus-driven malignant transformation, including viral-encoded miRs. Recently, a growing number of target genes regulated by EBV-miRs have been identified, ranging exemplarily from immunological relevant targets, like IL-12 by ebv-miR-BART1, ebv-miR-BART2, and ebv-miR-BHRF1-2,[65] TAP2 by ebv-miR-BART17 and ebv-miR-BHRF1-3,[27] MICB by ebv-miR-BART2-5p[28] and IL-1 receptor 1 (IL1R1) by ebv-miR-BHRF1-2-5p,[66] CXCL11 by ebv-miR-BHRF1-3,[67] to the pro-apoptotic genes PTEN[24] by ebv-miR-BART7-3p, PUMA by ebv-miR-BART20-5p and ebv-miR-BART8,[68] DAB2 by ebv-miR-BART1-3p,[69] DICE1 by ebv-miR-BART3[70] and TP53 by ebv-miR-BART-5-3p.[71]

In the process of malignant transformation, an altered glucose metabolism known as Warburg effect is considered as an early event. It represents the metabolic basis for tumor cell proliferation and growth and also contributes to immune evasion due to the local acidic pH of the tumor microenvironment.[72] EBV infection has been shown to reorganize B cell metabolic pathways to promote immortalization as a precursor to the development of lymphomas. An increased aerobic glycolysis in EBV-infected B cells was associated with an upregulation of lactate transporters.[73] This is in line with a statistically significant higher expression in several key enzymes involved into the early steps of glycolysis in BL cells when compared to the healthy B cells. These findings might encourage further studies to investigate the impact of EBV-miRs towards metabolic alterations in the process of malignant transformation.

Next to this, other immune-relevant putative EBV target genes were identified, which mainly represent genes involved in the viral pattern recognition, in the recruitment of immune effector cells (chemokines) towards virus-infected host cells and in anti-viral interferon (IFN) signaling. These mechanisms might be beneficial for early or latent stages of EBV infection. Furthermore, immune modulatory molecules were detected to be statistically significant enhanced upon EBV infection, including PD-L1 and CD155, both involved in the inhibitory signaling to immune effector cells thereby possibly contributing to immune evasion of EBV-infected cells. This is in line with reports demonstrating an induction of PD-L1 upon EBV infection with an involvement of LMP1.[74,75]

Another hallmark of cancer is the activated invasiveness and ability of metastasis.[51] The in silico part of this study identified statistically significant reduced gene expression of epithelial marker genes as well as the enhanced expression of mesenchymal marker genes in the BL cells when compared to the healthy B cells offering putative targets of the EBV-driven malignant transformation. However, these results might be also translatable to solid EBV-associated tumor entities, for which EMT and MET processes are of higher relevance than for BLs. In fact, various EBV-encoded molecules are known to enhance EMT, in particular ebv-miR-BART7-3p,[24] ebv-miR-BART8-3p,[76] ebv-miR-BART11,[77] ebv-miR-BART1[31] and LMP1.[78]

The combination of several other cancer hallmarks is reflected by the analysis of the gene expression pattern of factors involved in cell proliferation, cell cycle regulation as well as of known oncogenes. Many proliferation markers and known oncogenes exerted a statistically significant altered gene expression. A connection between EBV infection and higher NRAS gene expression has been published,[79] while an increased PDGFB expression has been reported in NPC, independently of EBV.[80] Together with the EBV-encoded LMP2A, SYK plays a role in cell migration of epithelial cells.[81] The proliferation marker and cell cycle regulator E2F1 is induced by the EBV-encoded BZLF1 protein.[82]

In this study, a total of 161 genes with relevance for immune surveillance/immune evasion as well as neoplastic transformation have been screened for their deregulation upon EBV infection. A total of 24 of these 161 genes analyzed have already been linked to EBV infection and for some cases, the underlying molecular mechanisms have been identified. Multifaceted cancer-related processes were detected to be deregulated upon EBV infection and the induced gene expression could be explained by regulatory mechanisms, such as enhanced PIP3-AKT and ERK signaling, caused by EBV-encoded molecules.

A follow-up in silico analysis also highlighted interesting EBV-miR candidates caused a statistically significant downregulation in 5 out of the 24 deregulated target genes, thereby representing putative direct targets of EBV-miRs. Interestingly, 3/5 downregulated genes are involved in IFN signaling. However, a direct interaction between these candidate target genes and the respective EBV-miRs has not yet been proven, while other EBV-mediated effects could directly or indirectly cause the observed downregulation. In addition, 41 tumor biological relevant target genes exerted an upregulation upon EBV infection with a majority of them already described in the literature in the context of EBV as exemplarily described above.

Next to this, it has been demonstrated that various EBV-encoded miRs are expressed and detectable at high levels during latency with heterogeneous levels and a tissue-specific expression patterns in infectious mononucleosis, B cell lymphomas, NPCs as well as in EBV-positive B cell lines, but not in healthy tumor-free tissues and EBV-negative B cell lines. Since these EBV-miRs were stable in FFPE specimen, they might represent suitable prognostic markers, in particular when other EBV detection methods lead to unclear results as e.g. demonstrated for the LMP1 staining associated with the phase of EBV infection (Table 2). It was an initial intention of this study to underline the possible usage of EBV-miRs as diagnostic and/or prognostic markers. A value of miRs as biomarkers has already been shown for several human endogenous miRs, such as e.g. the oncogenic miR-21.[83] However, there is still an urgent need to extend the knowledge to viral-encoded miRs.

The qPCR-based determination of the most abundant EBV-miRs might also represent a possible alternative for EBV detection and give in depth knowledge of EBV-associated diseases. Furthermore, the combination of EBV-miR expression analyses and in silico prediction studies appears to be very promising tool leading to the identification of novel EBV-miR target genes, while avoiding the actual limitations of online target prediction databases lacking EBV-miRs.

The selection of 161 analyzed genes led to the identification of 24 EBV-driven deregulated genes with a heterogenous chromosomal distribution, but with chromosome 12, 9, and 10 harboring the most deregulated genes. This is noteworthy, since abnormalities of these chromosomes frequently occur in hematopoietic malignancies, such as acute myeloid leukemia (AML). Translocations leading to the well-known oncofusion proteins BCR-ABL1 (t(9;22)), DEK-CAN (t(6;9)), and MN1-TEL (t(12;22)),[84] but also abnormalities in 11q23 have been linked to an unfavorable outcome in AML patients.[85]

Therapeutic strategies involving the targeting/inhibition of viral or endogenous oncogenic miRs have already reached preclinical animal studies and first clinical trials with promising results. These include the oncogenic miR-21 in melanoma, which was blocked by locked-nucleic acid miR inhibitors.[86] In the case of Miravirsen (miR-34a) and Cobomarsen (miR-155), miR-based therapies even entered clinical phase II and phase I trials for anti-tumoral therapies. In addition, an increasing number of patents are registered based on miRs as therapeutics against several viruses including the Ebola virus, Hepatitis C virus and the human enterovirus 71. Furthermore, miR-based therapies against cardiovascular diseases, renal diseases and other diseases have entered clinical studies, which have been recently reviewed by Chakraborty and co-authors.[87]

Conclusions

This approach employs an in silico analysis to screen for novel so far unidentified EBV-targeted tumor biology relevant genes in the host cell by circumvention of the limited possibilities of the available databases in the context of EBV-encoded miRs. Numerous EBV-miRs were analyzed in different EBV-positive human (tumor) tissues/cells, leading to the identification of certain abundantly expressed EBV-miRs suitable for EBV detection, which might also have prognostic potential with an advantage to the conventional methods only analyzing one marker. Despite a limited number of genes were used for screening, five novel so far unidentified putative target genes for EBV-encoded miRs were identified, whose gene products were involved in tumor immunologic and tumor biologic related processes. An increased knowledge of the crucial role of EBV-miRs in EBV-associated malignancies and their targets does not only provide insights into the underlying mechanisms, but also represents a rational for the development of therapeutic strategies for these diseases.