Comprehensive profiling of functional Epstein-Barr virus miRNA expression in human cell lines
BMC Genomics volume 17, Article number: 644 (2016)
Epstein-Barr virus (EBV) establishes lifelong infections in its human host. The virus is associated with a broad range of malignancies of lymphoid and epithelial origin, including Burkitt’s lymphoma, post-transplant lymphoproliferative disease, nasopharyngeal carcinoma and gastric carcinoma. During the latent phase of its life cycle, EBV expresses more than 40 mature miRNAs that are highly abundant in tumor cells and may contribute to oncogenesis. Although multiple studies have assessed the relative expression profiles of EBV miRNAs in tumor cells, data linking these expression levels to functional target knockdown are mostly lacking. Therefore we set out to systematically assess the EBV miRNA expression levels in EBV+ tumor cell lines, and correlate this to their functional silencing capacity in these cells.
We provide comprehensive EBV miRNA expression profiles of the EBV+ cell lines C666-1 (nasopharyngeal carcinoma), SNU-719 (gastric carcinoma), Jijoye (Burkitt’s lymphoma), and AKBM (Burkitt’s lymphoma) and of EBV− cells ectopically expressing the BART miRNA cluster. By deep sequencing the small RNA population and conducting miRNA-reporter experiments to assay miRNA potency, we were able to compare the expression profiles of the EBV miRNAs with their functional silencing efficacy. We observe a strong correlation between miRNA expression levels and functional miRNA activity. There is large variation in expression levels between EBV miRNAs in a given cell line, whereas the relative expression profiles are well maintained between cell lines. Furthermore, we show that miRNA arm selection bias is less pronounced for gamma-herpesvirus miRNAs than for human miRNAs.
We provide an in depth assessment of the expression levels and silencing activity of all EBV miRNAs in B- and epithelial cell lines of different latency stages. Our data show a good correlation between relative EBV miRNA expression levels and silencing capacity, and suggest preferential processing of particular EBV miRNAs irrespective of cell-type. In addition to encoding the largest number of precursor miRNAs of all human herpesviruses, EBV expresses many miRNAs precursors that yield two functional miRNA strands, rather than one guide strand and a non-functional passenger strand. This reduced strand bias may increase the size of the EBV miRNA targetome.
MicroRNAs (miRNAs) are small RNA molecules of ~22 nucleotides in length that regulate gene expression by binding to target mRNAs in RNA-induced silencing complexes (RISC). Mechanisms for miRNA-mediated target gene downregulation include enhanced mRNA degradation, suppression of translation and in rare cases slicing of the mRNA at the target site (reviewed in ). In the last decade, many virus-encoded miRNAs have been identified, predominantly within the genomes of large dsDNA viruses such as herpesviruses. The two members of the gamma-herpesvirus subfamily that infect humans, Kaposi Sarcoma Associated Herpesvirus (KSHV, HHV-8) and Epstein-Barr virus (EBV, HHV-4), express large clusters of miRNAs during their latent phase (reviewed in ). Since very few viral proteins are expressed during latency, the expression of these miRNAs greatly expands the number of gene products available to manipulate the host cell. Indeed, EBV miRNAs play important roles in regulating EBV-induced B cell transformation [3–5] and downmodulate pro-apoptotic proteins such as PUMA , Bim , Bid , caspase 3  and others . Additionally, the EBV miRNAs contribute to immune evasion by targeting the NK cell ligand MICB , the chemokine CXCL11  and the inflammasome component NLRP3 .
EBV establishes a lifelong infection in over 90 % of the adult human population, where it predominantly resides in a latent stage in the memory B cell compartment. During its life cycle, the virus passes through different phases that are characterized by distinct gene expression patterns . After primary B cell infection, EBV gradually shuts down its protein expression. Initially, during the latency III stage, all EBNA (EBV nuclear antigens) and LMP (latent membrane protein) proteins are expressed. This is followed by the latency II stage where only EBNA1 and LMP proteins are expressed. During the subsequent latency I stage, protein expression is restricted to EBNA-1 alone. Through occasional reactivation and possibly lytic replication in epithelial cells of the oropharynx, EBV can spread to new hosts . Although primary EBV infection usually proceeds asymptomatically, it is the causative agent for infectious mononucleosis. Furthermore, EBV infection is associated with a number of malignancies of epithelial, B cell and (occasionally) NK/T cell origin . EBV-associated tumors display one of the three gene expression profiles described above. Burkitt’s lymphoma (BL) for instance, is characterized by a latency I program, although EBV frequently switches to latency III in vitro . Hodgkin’s lymphoma (HL) and nasopharyngeal carcinoma (NPC) follow the latency II program, whereas gastric carcinoma (GC) is characterized by a latency I expression pattern, although occasionally LMP2A can be detected as well [16–18]. In EBV-associated lymphomas of immunosuppressed individuals, EBV expresses the largest set of latent genes (latency III) . Also B lymphocytes transformed with EBV in vitro (lymphoblastoid cell lines, LCLs) display the latency III gene expression pattern .
Intriguingly, whereas viral protein expression is strictly regulated during the EBV latency stages, viral miRNAs are expressed throughout all stages of latency. EBV miRNAs are expressed in clusters and are named BART and BHRF1 miRNAs [21–23]. The BART miRNAs are located in introns of the latent BART transcripts (BARTs) that are expressed during all latency phases. Upon induction of the lytic cycle, the expression of the BARTs is increased, which coincides with a modest increase in mature BART miRNA expression levels . The BHRF1 miRNAs lie upstream (miR-BHRF1-1) or downstream (miR-BHRF1-2, miR-BHRF1-3) the protein-coding BHRF1 gene and are mainly processed from latent EBNA transcripts that are highly expressed in latency III B cells. EBV reactivation results in expression of lytic BHRF1 mRNA from which miR-BHRF1-2 and miR-BHRF1-3 can be processed [24, 25]. The relative abundance of EBV miRNAs has been assessed in several human tissues and cell lines via deep sequencing, microarray, or qPCR techniques [8, 9, 24, 26–32]. Although insightful, these evaluations are limited to the relative expression profiles of EBV miRNAs without correlating this to their biological silencing potential. EBV miRNA targets have been identified by various approaches, including PAR-CLIP studies [9, 32], but a comprehensive analysis of the activities of all EBV miRNAs is lacking. Here, we set out to connect these two parameters experimentally. For this, we monitored miRNA activity of both strands from each EBV-encoded pre-miRNA in four different cell lines that cover the main latency stages of EBV. In parallel we performed small RNA deep sequencing in these same lines and correlated their relative expression level to functional target knockdown. Our data show a strong correlation between miRNA expression level and functional target knockdown. EBV miRNA expression profiles are highly consistent between cell lines, where miRNAs abundantly expressed in one cell line are also abundantly expressed in other EBV+ cell lines. Expression of the BART cluster in EBV-negative cells maintains the same hierarchy. Intriguingly, we noticed that there is considerably less strand bias in the expression of viral miRNAs as compared to human miRNAs: viral miRNAs have a higher incidence of expressing both miRNA strands in cells, thereby potentially enhancing the regulatory capacity of these miRNA genes.
Results and discussion
Functional expression profiling using miRNA reporters
In order to study the functional expression of EBV miRNAs in infected cells, lentiviral reporter constructs containing a single perfect miRNA binding site downstream of the fluorescent mCherry gene (Fig. 1a) were designed for all mature EBV miRNAs based on their sequences annotated in miRBase . Although many miRNA hairpins are known to favor selection of one arm of the duplex over the other (the guide strand versus the star strand), both arms of a miRNA may be functional and capable of downregulating target mRNAs. We aimed at generating a comprehensive profile of the functional expression of all mature EBV miRNAs including the ones lacking miRBase annotation. To assess the functional expression of the latter ones, we designed miRNA reporters based on previous deep sequencing studies  or we introduced the complementary sequence of the entire miRNA arm including several adjacent nucleotides as a potential target site in our lentiviral mCherry reporter constructs. A complete list of the EBV miRNA sensor sequences is presented in Additional file 1.
We selected four EBV+ cell lines for functional screening: SNU-719 (type 1 EBV, gastric carcinoma, latency I/II, ), C666-1 (type 1 EBV, nasopharyngeal carcinoma, latency II), Jijoye (type 2 EBV, Burkitt’s lymphoma, latency III) and AKBM (type 1 EBV, Burkitt’s lymphoma, latency I, Akata strain). Additionally, in order to study EBV miRNA expression outside the context of virus infection, we analyzed EBV-negative HK1 (nasopharyngeal epithelial cell line) ectopically expressing the complete BART miRNA cluster but lacking BART2 (“HK1-miR-BART”). We next measured the functional downregulation of each miRNA reporter in these five cell lines (Fig. 1b–f). As expected, we observed differential downregulation of the miRNA reporter constructs, and these assessments were highly reproducible between experiments (Additional file 2). The majority of the BART miRNAs induced potent downregulation of reporters in the carcinoma cell lines (Fig. 1b–c). These findings are consistent with high expression levels of BART miRNAs in epithelial cells reported previously [8, 30, 36]. BART miRNA-mediated sensor downregulation in HK1-miR-BART cells was similar to that in EBV+ epithelial cells (Fig. 1f). For all miRNAs, at least one of the two strands induced potent target knockdown (>70 % repressed in the carcinoma cell lines), with the exception of miR-BART21 and miR-BART20 for which both strands induced only minor target downregulation (±50 %). In agreement with our data, large variations in expression levels between BART miRNAs have been observed previously [21, 29, 37, 38]; the mechanism behind this regulation remains to be determined.
The BHRF1 miRNAs are processed from EBNA-encoding transcripts initiated at the Cp or Wp promoter [21, 24], which are active during latency III. In addition, miR-BHRF1-2 and 3 can be processed from the 3’UTR of the BHRF1 mRNA, a lytic transcript that has also been detected in latent cells with an active Wp promoter . Consistent with this, we barely see functional expression of BHRF1 miRNAs in carcinoma cell lines C666-1 and SNU-719, in which the Qp promoter drives expression of EBNA1 and the Cp and Wp promoter are not active (Fig. 1b–c). In contrast, in the latency III Jijoye cells, miR-BHRF1-2-5p, miR-BHRF1-2-3p and miR-BHRF1-3-5p potently downregulated their respective reporters (Fig. 1d). The miR-BHRF1-1-3p and miR-BHRF1-1-5p miRNAs do not show silencing activity in Jijoye cells. This is likely the result of a single nucleotide change in the 3p arm of the miR-BHRF1-1 pre-miRNA , which may preclude proper miR-BHRF1-1 expression.
The reporters for most of the miRNA strands absent from miRBase (miR-BART12-5p, miR-BART15-5p, miR-BART22-5p, miR-BHRF1-1-3p and miR-BHRF1-3-3p) were not or barely downregulated in the various cell lines tested, indicating that the registration of functional mature EBV miRNAs in this database is quite complete. However, although the miR-BART16-3p reporter was clearly downregulated in most cell lines, this miRNA is not annotated in miRBase and has not been analyzed in any of the previously published qPCR profiling studies of EBV+ samples.
A strong correlation in BART miRNA reporter downregulation was observed between SNU-719, C666-1 and Jijoye (Fig. 2), i.e. in all three cell lines the same hierarchy of miRNA potency was observed between individual mature miRNAs. Remarkably, this relative expression profile was also maintained in the EBV-negative HK1 cell line expressing the BART miRNAs. This indicates that the main determinants for BART pre-miRNA processing, miRNA stability, RISC incorporation and target binding lie in the primary miRNA sequence (pri-miRNA) and are not affected significantly by viral sequences outside of the miRNA cluster or by cell-type specific expression of host cell genes. These findings are consistent with a study by Qiu et al.,  in which EBV miRNA expression patterns were highly similar in tumor biopsies from various origins. Although we observed comparable miRNA expression patterns between cell lines, overall the BART miRNA reporters were less potently downregulated in the Jijoye cells as compared to the epithelial cell lines. This observation corresponds with a study by Cai et al. , where BART transcript levels were reduced in Jijoye cells as compared to C666-1 cells.
AKBM cells displayed a much-reduced potency of EBV miRNA activity as compared to the other cell lines; the majority of mature miRNAs induced only minor sensor downregulation (ranging from 0 to 25 %), whereas only miR-BART2-3p and miR-BART3-5p reduced mCherry expression by ±50 % (Fig. 1e, Additional file 3). Since these reporters acted distinctively as compared to the majority of EBV miRNA sensors, we tested whether they could be targeted by human instead of viral miRNAs. Although miR-BART2-3p and miR-BART3-5p reporters were not downregulated in EBV-negative HK1 cells (data not shown), we did observe a marked downregulation of the BART3-5p sensor in EBV-negative subclones of the Akata line which have lost the EBV genome (Additional file 4). This suggests that the observed regulation of this reporter in EBV-positive AKBM cells is likely caused by human miRNAs expressed in these cells, and not by viral miRNAs. The miR-BART2-3p reporter, however, was not downregulated in EBV-negative Akata cells, indicating that this miRNA may be functional in AKBM cells. It is unclear why miR-BART2-3p is more potent in AKBM cells as compared to other EBV miRNAs expressed in these cells.
Expression profiling using deep sequencing
To allow correlation of our functional EBV miRNA profiling to mature miRNA expression levels, we next assessed the relative abundance of all EBV miRNAs by deep sequencing. This yielded 32–41 million sequencing reads per sample, after filtering of adapter sequences and low-quality reads (Table 1, Additional file 5). Reads were aligned to the human and EBV genome. The average read length was ~22 nucleotides, of which the majority mapped to miRNA sequences (75–97 %). The proportion of EBV reads ranged from 23.5 % in Jijoye to only 1.4 % in AKBM (Fig. 3a), suggesting that the poor miRNA reporter downregulation in the Akata-derived AKBM cells reflects low expression levels of the given miRNAs. To ensure that the relatively low expression of EBV miRNAs was not caused by loss of EBV from a large percentage of cells, we constructed a reporter vector in which 5 EBV miRNAs target sites were cloned downstream the mCherry reporter gene. Upon introduction of this ‘super-sensor’ in AKBM cells, we observed silencing in the majority of the cells which was not seen in EBV-negative subclones of the Akata cell line 2A8 and AK31 (Additional file 6), showing that relatively low EBV miRNA expression was not caused by absence of EBV from the AKBM cells.
For the C666-1, SNU-719, and Jijoye cells, there was a clear inverse correlation between relative read counts (reads per million reads, RPM) and miRNA reporter expression (Fig. 3b). However, the correlation between these readouts within a single cell line (Fig. 3b) was less pronounced than the correlation of either read counts (Additional file 7) or reporter downregulation between two cell lines, e.g. the C666-1 and SNU-719 cells (Fig. 2). This difference may be explained by preferential deep sequencing of certain miRNA sequences; studies have indicated that miRNA sequencing bias can be introduced in the adapter ligation procedure, thereby favoring certain mature miRNA sequences over others . Moreover, it is unlikely that the potency of individual miRNAs depends solely on the relative abundance of this miRNA. Indeed, the relative location of the miRNA target site in the 3’UTR impacts the potency of regulation . In our experiments however, the distance between the mCherry reporter stop codon and the single miRNA target site is identical for all reporters. Besides target-site location, also seed-pairing stability (strength of the interaction) and target site abundance have been shown to influence miRNA activity .
We observed less miRNA reporter downregulation in Jijoye cells as compared to epithelial cell lines, although sequence read counts for given miRNAs were comparable (Figs. 1, and 3b, and Additional file 5). However, since we only assessed the relative EBV miRNA abundance (reads per million), it is possible that the Jijoye cells may actually have low absolute levels of EBV miRNAs. Indeed, previous studies report lower absolute BART miRNA expression levels in Jijoye as compared to C666-1 [21, 37]. Similarly, we observed roughly 2–3 fold lower relative EBV miRNA reads in the HK1-miR-BART cells as compared to the C666-1 and SNU-719 cells (Fig. 3a), although the potency of downregulation of the miRNA sensors was comparable between these three lines (Fig. 1). It is possible that the absolute EBV miRNA levels is similar between these lines.
Deep sequencing studies revealed that miRNAs can display extensive length and sequence variation . These variants, termed isomiRs, are generated by differential Drosha/Dicer processing, trimming by exonucleases, non-templated nucleotide additions, and RNA editing events . Indeed, we observed a large variety of isomiRs for all EBV miRNAs in each of the EBV+ cell lines analyzed. For many miRNAs less than 50 % of the sequenced reads had the exact canonical sequence as reported in miRBase (Additional file 5). Additionally, the most abundant sequence for several EBV miRNAs did not match the annotated sequence deposited in miRBase. The majority of isomiRs displayed 3’end variation, which is not expected to severely impact miRNA activity due to preservation of the 5’ end seed sequence. Nevertheless, the most abundant miR-BART10-3p sequence lacks a single 5’ nucleotide as compared to the canonical sequence in the EBV+ cell lines. Besides for miR-BART10-3p, high frequencies of isomiRs with alternative 5’ ends were also apparent for miR-BART2-3p, miR-BART16-3p and miR-BART5-5p (Additional file 5). Of these, only miR-BART10-3p and miR-BART5-5p represent the dominant arm of the miRNA duplex. For miR-BART5-5p, 38–60 % of the reads have a single nucleotide deletion at their 5’ end as compared to miRBase (Additional file 5). It is likely that these 5’ end variants will have an altered gene target repertoire as compared to the canonical miRNA sequences. Detection of the 5’ end variability prompted us to assess whether our miRNA reporters could be targeted by these isomiRs. Fortunately, only a minority of the reads had an elongated 5’ end that could obstruct activity towards our reporter construct. The relatively abundant 5’ variants of miR-BART10-3p, miR-BART2-3p, and miR-BART5-5p were truncated at their 5’ termini, rendering the reporters used in these studies still functional. For miR-BART16-3p, however, ~40 % of the sequences possessed a 2 nt extension at the 5’ end for which no complementary bases were present in our miRNA reporter construct. It is therefore possible that we underestimate the functional activity of this specific miRNA.
Non-templated nucleotides can be added to precursor or mature miRNA 3’ termini by terminal nucleotidyl transferases (TNTs). The functional implications of these modifications are not fully understood, although miRNA uridylation and adenylation have been reported to impact miRNA stability [46–49] and potency [50, 51], and may act as a molecular signal to facilitate miRNA transport within and between cells . We observed low levels of non-templated nucleotide additions (NTAs) to EBV miRNAs in our cell lines (ranging from 0 to 24 %) where adenine and uracil additions were most frequent. Intriguingly, we observed a clear correlation between the frequency of NTAs for given miRNAs between cell lines (Fig. 4a and b), where some EBV miRNAs were clearly preferentially adenylated or uridylated. In line with our findings, previous reports describe sequence specificity of terminal nucleotidyl transferases [48, 51, 53]. Whether these NTAs are important for EBV miRNA function remains to be determined.
Reduced strand bias of gamma-herpesvirus miRNAs
Upon loading of a miRNA duplex onto Argonaute, one miRNA strand/arm (the guide RNA) is selected whereas the other one (passenger strand) is degraded . For many miRNAs, there is a strong bias towards selection of one strand over the other, which results in profound differences in abundance of both miRNA strands within cells . Although we saw a clear strand bias for human miRNAs, this phenomenon appeared much weaker for EBV miRNAs. For approximately half of the BART miRNAs, one of both arms was clearly dominant, whereas in a third of cases both strands potently downregulated their reporter (Fig. 1). This same feature was apparent in the deep-sequencing data: where for human miRNAs the median ratio between the most and least abundant miRNA strand was 1:186, for EBV miRNAs this was only 1:11 (Fig. 5a, left two bars). For a quarter of the EBV miRNAs there was a less than twofold difference between both strands, which applied to only 8 % of the human miRNAs. This difference in miRNA arm bias between human and EBV miRNAs was apparent for both high and low expressed miRNAs, as subdividing the miRNAs based on abundance yielded the same trend (Fig. 5a, right six bars). This phenomenon was not limited to C666-1 cells as we observed the same in SNU-719, Jijoye, and AKBM cells in our dataset and for Akata-LCL and BC-1 cells in sequencing data from other research groups (Fig. 5b).
We next extended our strand-bias analysis to non-EBV gamma-herpesvirus miRNAs (Fig. 5c). Intriguingly, also KSHV miRNAs expressed in BC-1 cells displayed reduced strand bias as compared to host miRNAs and the same was apparent for rhesus lymphocryptovirus (rhLCV, an EBV ortholog) miRNAs in rhesus monkey LCLs (Fig. 5b). Hence, it appears that gamma-herpesvirus miRNAs have retained or acquired the ability of both strands to functionally downregulate target genes.
The main feature determining arm selection is the thermodynamic stability of the miRNA duplex, whereby the miRNA arm with the least stable 5’ end is typically the guide strand that is loaded into RISC [56–58]. Also the secondary structure of the miRNA  and the identity of its first nucleotide [58, 60, 61] have been described to impact strand selection. How this reduced strand bias for gamma-herpesvirus miRNAs is obtained or retained and why this occurs remains to be determined. In principle, reducing strand bias may allow EBV to expand the number of target genes it regulates. Increasing the expression of miRNAs from one to two strands can be viewed as a genomic space-saving approach to enlarge regulatory capacity. This may be especially favorable for viruses since they are more constrained in their genome size than their host.
We have provided comprehensive functional and relative expression profiles for all EBV miRNAs in four EBV+ cell lines representing different latency stages. In general, there is good concordance between the relative expression levels of miRNAs and their potential to silence target genes. The dominance of particular miRNAs over others is shared by all four cell lines and an EBV− cell line in which the BART miRNAs were introduced through lentiviral transduction. This observation argues against an important role for cell-specific or virus-specific factors in regulating miRNA activity or abundance (apart from promoter usage), but rather suggests that the (pri-) miRNA sequence determines functional expression levels. Our deep sequencing datasets provide in depth information on EBV miRNA sequence diversity, including variation that affects miRNA target specificity by 5’end heterogeneity. Interestingly, EBV not only expresses the highest numbers of miRNAs of all human viruses; it also has increased regulatory capacity by functionally expressing both arms for a significant number of miRNAs.
For the cloning of miRNA reporters, single perfect miRNA binding sites were introduced 29 bp downstream the mCherry reporter gene in the pSicoR-EF1a-PuroR-T2A-mCherry vector (described elsewhere, manuscript in preparation). See Additional file 1 for sequences cloned. Of note: the deep-sequencing analysis identified a single nucleotide change (position 17, T changed to C) in the mature miR-BART19-5p miRNA sequence in the C666-1 strain. Therefore, the miR-BART19-5p miRNA sensor used in the C666-1 cells to monitor miRNA activity does not fully match the miR-BART19-5p sequence expressed in these cells.
We ectopically expressed EBV miRNAs in EBV− HK1 cells from the pLenti-Blast-eGFPintron vector (described elsewhere, manuscript in preparation). The miRNAs were cloned into an intron located in an eGFP marker gene driven by the EF1a promoter. Briefly, the EBV BART region was PCR amplified from genomic EBV DNA in three different PCR reactions (BART3 to 6, 21 to 18 and 7 to 14), which were subsequently combined into the pLenti-Blast-eGFPintron vector. All constructs were verified by Sanger sequencing.
Cells and lentiviral transduction
Jijoye and 293 T human embryonic kidney cells were obtained from the ATCC (American Type Culture Collection). HK1 cells  were kindly provided by Prof. Dr. G. Tsao (The University of Hong Kong, Hong Kong). SNU-719 cells (originally derived from the Korean Cell Line Bank) and C666-1 were a kind gift from Prof. Dr. J. Middeldorp (VU University Medical Center, Amsterdam). As C666-1 cells can lose LMP1 and LMP2A expression in culture, the LMP1/2A expression status was verified by TaqMan qPCR which showed expression of both viral transcripts in our cells (data not shown). The EBV+ AKBM cell line  and EBV− AK31 and 2A8 cell lines are derived from the Japanese Burkitt’s lymphoma cell line Akata. Cell lines were maintained in RPMI (Roswell Park Memorial Institute) 1640 supplemented with heat-inactivated 10 % fetal calf serum, 2 mM glutamine and antibiotics. Tissue culture flasks for culturing C666-1 cells were coated with fibronectin (Sigma F1141, ). Lentiviruses were produced via standard lentiviral production protocols using third-generation packaging vectors. Transductions of EBV+ cell lines were performed by spin infections at 1000 g for 1.5 h at 33 °C in the presence of 4 ug/ml polybrene and lentivirus. To obtain a pure population of BART miRNA expressing HK1 cells, lentivirally transduced cells were FACS sorted on a BD FACS Aria II based on eGFP expression.
MiRNA reporter assay
Cells were lentivirally transduced at low MOI to achieve single lentiviral integrations, allowed to recover for ± 1 week, fixed, and analyzed by flow cytometry on a BD Canto II flow cytometer. FlowJo software (Treestar) was used for analyzing the flow cytometric data. The gating strategy consisted of selecting live cells (based on forward and side scatter) and then gating on mCherry-positive cells (thereby gating out untransduced cells). To calculate relative reporter expression, the mCherry geometric means of cells harbouring the miRNA reporters were normalized against those of cells harbouring the empty control vector.
Construction of small RNA libraries and deep sequencing
RNA was isolated using Trizol according to manufacturer’s instructions (Life Technologies). Integrity of the RNA was verified using the Agilent 2100 Bioanalyzer (RIN scores were 10.0). Sequencing of small RNA was performed by BGI (Hong Kong, China), briefly: 3’ and 5’ adapters were ligated to size-selected RNA molecules (~18–30 nt length), followed by reverse transcription using SuperScript II (Invitrogen). cDNA was PCR amplified, after which clusters were generated using an Illumina cBot cluster amplification system. Single-end sequencing by synthesis was performed on the HiSeq 2000 platform using the TruSeq SBS Kit-HS V3. Good quality sequence reads were trimmed to remove adapter sequences and analyzed by using the sRNAbench 0.9 package . Reads were aligned to the human (reference genome GRCh37) and EBV genome (NC_007605). Initially, no mismatches were allowed between the Illumina reads and the reference sequences. However, as we identified a single nucleotide change (position 17, T changed to C) in the mature BART19-5p miRNA sequence in the C666-1 EBV strain (as reported before by Chen et al. ), read count numbers for this miRNA were adjusted to reflect the correct miR-BART19-5p reads. We did not detect additional polymorphisms in any of the analyzed cell lines. Correlation coefficients were calculated using Graphpad Prism.
MiRNA arm ratio calculations
miRNA arm ratios were calculated by dividing the number of reads of the most abundant arm by that of the least abundant arm, to obtain numbers between 1 (1 : 1 ratio) and infinite (if one arm was undetectable). Data from published studies were obtained from the NCBI Sequence Read Archive (http://www.ncbi.nlm.nih.gov/sra): Akata-LCL: SRR1003028 , BC-1: SRR343332 , rLCV-infected 211-98 rhesus macaque rLCL: run SRR1003031 and run SRR1003028, rLCV-infected 309-98 rhesus macaque rLCL: run SRR1003035 and run SRR1003033 . Deep sequencing reads were analyzed by using the sRNAbench 0.9 package, single nucleotide mismatches to reference sequences were allowed in the analysis.
BL, Burkitt’s lymphoma; EBV, Epstein-Barr virus; GC, gastric carcinoma; HL, Hodgkin’s lymphoma; KSHV, Kaposi’s sarcoma associated herpesvirus; LCL, lymphoblastoid cell line; mRNA, messenger RNA; miRNA, microRNA; NPC, nasopharyngeal carcinoma; rhLCV, rhesus lymphocryptovirus; RISC, RNA induced silencing complex.
Ameres SL, Zamore PD. Diversifying microRNA sequence and function. Nat Rev Mol Cell Biol. 2013;14(8):475–88.
Zhu Y, Haecker I, Yang Y, Gao SJ, Renne R. gamma-Herpesvirus-encoded miRNAs and their roles in viral biology and pathogenesis. Curr Opin Virol. 2013;3(3):266–75.
Feederle R, Linnstaedt SD, Bannert H, Lips H, Bencun M, Cullen BR, et al. A viral microRNA cluster strongly potentiates the transforming properties of a human herpesvirus. PLoS Path. 2011;7(2), e1001294.
Seto E, Moosmann A, Gromminger S, Walz N, Grundhoff A, Hammerschmidt W. Micro RNAs of Epstein-Barr virus promote cell cycle progression and prevent apoptosis of primary human B cells. PLoS Pathog. 2010;6(8), e1001063.
Vereide DT, Seto E, Chiu YF, Hayes M, Tagawa T, Grundhoff A, et al. Epstein-Barr virus maintains lymphomas via its miRNAs. Oncogene. 2014;33(10):1258–64.
Choy EY, Siu KL, Kok KH, Lung RW, Tsang CM, To KF, et al. An Epstein-Barr virus-encoded microRNA targets PUMA to promote host cell survival. J Exp Med. 2008;205(11):2551–60.
Marquitz AR, Mathur A, Nam CS, Raab-Traub N. The Epstein-Barr Virus BART microRNAs target the pro-apoptotic protein Bim. Virology. 2011;412(2):392–400.
Shinozaki-Ushiku A, Kunita A, Isogai M, Hibiya T, Ushiku T, Takada K, Fukayama M: Profiling of Virus-encoded MicroRNAs in Epstein-Barr Virus-associated Gastric Carcinoma and their Roles in Gastric Carcinogenesis. J Virol. 2015;89(10):5581-91.
Kang D, Skalsky RL, Cullen BR. EBV BART MicroRNAs Target Multiple Pro-apoptotic Cellular Genes to Promote Epithelial Cell Survival. PLoS Pathog. 2015;11(6), e1004979.
Nachmani D, Stern-Ginossar N, Sarid R, Mandelboim O. Diverse herpesvirus microRNAs target the stress-induced immune ligand MICB to escape recognition by natural killer cells. Cell Host Microbe. 2009;5(4):376–85.
Xia T, O’Hara A, Araujo I, Barreto J, Carvalho E, Sapucaia JB, et al. EBV microRNAs in primary lymphomas and targeting of CXCL-11 by ebv-mir-BHRF1-3. Cancer Res. 2008;68(5):1436–42.
Haneklaus M, Gerlic M, Kurowska-Stolarska M, Rainey AA, Pich D, Mcinnes IB, et al. Cutting edge: miR-223 and EBV miR-BART15 regulate the NLRP3 inflammasome and IL-1beta production. J Immunol. 2012;189(8):3795–9.
Longnecker RM, Kieff E, Cohen JI. Epstein-Barr Virus. In: Knipe DM, Howley PM, editors. Field’s Virology. Philadelphia: Lippincott Williams & Wilkins; 2013. p. 1898–959.
Hadinoto V, Shapiro M, Sun CC, Thorley-Lawson DA. The dynamics of EBV shedding implicate a central role for epithelial cells in amplifying viral output. PLoS Pathog. 2009;5(7), e1000496.
Gregory CD, Rowe M, Rickinson AB. Different Epstein-Barr virus-B cell interactions in phenotypically distinct clones of a Burkitt’s lymphoma cell line. J Gen Virol. 1990;71(Pt 7):1481–95.
Sugiura M, Imai S, Tokunaga M, Koizumi S, Uchizawa M, Okamoto K, et al. Transcriptional analysis of Epstein-Barr virus gene expression in EBV-positive gastric carcinoma: unique viral latency in the tumour cells. Br J Cancer. 1996;74(4):625–31.
Pathmanathan R, Prasad U, Sadler R, Flynn K, Raab-Traub N. Clonal proliferations of cells infected with Epstein-Barr virus in preinvasive lesions related to nasopharyngeal carcinoma. N Engl J Med. 1995;333(11):693–8.
Pallesen G, Hamilton-Dutoit SJ, Rowe M, Young LS. Expression of Epstein-Barr virus latent gene products in tumour cells of Hodgkin’s disease. Lancet. 1991;337(8737):320–2.
Young L, Alfieri C, Hennessy K, Evans H, O’Hara C, Anderson KC, et al. Expression of Epstein-Barr virus transformation-associated genes in tissues of patients with EBV lymphoproliferative disease. N Engl J Med. 1989;321(16):1080–5.
Rickinson AB, Kieff E: Epstein-Barr Virus. In: Fields Virology. Edited by Knipe DM, Howley PM, vol. 5: Lippincott Williams & Wilkins; 2007: 2655–2700.
Cai X, Schafer A, Lu S, Bilello JP, Desrosiers RC, Edwards R, et al. Epstein-Barr virus microRNAs are evolutionarily conserved and differentially expressed. PLoS Pathog. 2006;2(3), e23.
Pfeffer S, Zavolan M, Grasser FA, Chien M, Russo JJ, Ju J, et al. Identification of virus-encoded microRNAs. Science. 2004;304(5671):734–6.
Grundhoff A, Sullivan CS. Virus-encoded microRNAs. Virology. 2011;411(2):325–43.
Amoroso R, Fitzsimmons L, Thomas WA, Kelly GL, Rowe M, Bell AI. Quantitative studies of Epstein-Barr virus-encoded microRNAs provide novel insights into their regulation. J Virol. 2011;85(2):996–1010.
Xing L, Kieff E. Epstein-Barr virus BHRF1 μ- and stable RNAs during latency III and after induction of replication. J Virol. 2007;81(18):9967–75.
Yang HJ, Huang TJ, Yang CF, Peng LX, Liu RY, Yang GD, et al. Comprehensive profiling of Epstein-Barr virus-encoded miRNA species associated with specific latency types in tumor cells. Virol J. 2013;10(1):314.
Zeng Z, Huang H, Huang L, Sun M, Yan Q, Song Y, et al. Regulation network and expression profiles of Epstein-Barr virus-encoded microRNAs and their potential target host genes in nasopharyngeal carcinomas. Sci China Life Sci. 2014;57(3):315–26.
Marquitz AR, Mathur A, Chugh PE, Dittmer DP, Raab-Traub N. Expression profile of microRNAs in Epstein-Barr virus-infected AGS gastric carcinoma cells. J Virol. 2014;88(2):1389–93.
Pratt ZL, Kuzembayeva M, Sengupta S, Sugden B. The microRNAs of Epstein-Barr Virus are expressed at dramatically differing levels among cell lines. Virology. 2009;386(2):387–97.
Cosmopoulos K, Pegtel M, Hawkins J, Moffett H, Novina C, Middeldorp J, et al. Comprehensive profiling of Epstein-Barr virus microRNAs in nasopharyngeal carcinoma. J Virol. 2009;83(5):2357–67.
Qiu J, Cosmopoulos K, Pegtel M, Hopmans E, Murray P, Middeldorp J, et al. A novel persistence associated EBV miRNA expression profile is disrupted in neoplasia. PLoS Pathog. 2011;7(8), e1002193.
Skalsky RL, Corcoran DL, Gottwein E, Frank CL, Kang D, Hafner M, et al. The viral and cellular microRNA targetome in lymphoblastoid cell lines. PLoS Pathog. 2012;8(1), e1002484.
Kozomara A, Griffiths-Jones S. miRBase: annotating high confidence microRNAs using deep sequencing data. Nucleic Acids Res. 2014;42(D1):D68–73.
Chen SJ, Chen GH, Chen YH, Liu CY, Chang KP, Chang YS, et al. Characterization of Epstein-Barr virus miRNAome in nasopharyngeal carcinoma by deep sequencing. PLoS One. 2010;5(9), e12745.
Oh ST, Seo JS, Moon UY, Kang KH, Shin DJ, Yoon SK, et al. A naturally derived gastric cancer cell line shows latency I Epstein-Barr virus infection closely resembling EBV-associated gastric cancer. Virology. 2004;320(2):330–6.
Kim Do N, Chae HS, Oh ST, Kang JH, Park CH, Park WS, et al. Expression of viral microRNAs in Epstein-Barr virus-associated gastric carcinoma. J Virol. 2007;81(2):1033–6.
Edwards RH, Marquitz AR, Raab-Traub N. Epstein-Barr virus BART microRNAs are produced from a large intron prior to splicing. J Virol. 2008;82(18):9094–106.
Grundhoff A, Sullivan CS, Ganem D. A combined computational and microarray-based approach identifies novel microRNAs encoded by human gamma-herpesviruses. RNA. 2006;12(5):733–50.
Kelly GL, Long HM, Stylianou J, Thomas WA, Leese A, Bell AI, et al. An Epstein-Barr virus anti-apoptotic protein constitutively expressed in transformed cells and implicated in burkitt lymphomagenesis: the Wp/BHRF1 link. PLoS Pathog. 2009;5(3), e1000341.
Cai X, Cullen BR. Transcriptional origin of Kaposi’s sarcoma-associated herpesvirus microRNAs. J Virol. 2006;80(5):2234–42.
Hafner M, Renwick N, Brown M, Mihailovic A, Holoch D, Lin C, et al. RNA-ligase-dependent biases in miRNA representation in deep-sequenced small RNA cDNA libraries. RNA. 2011;17(9):1697–712.
Grimson A, Farh KK, Johnston WK, Garrett-Engele P, Lim LP, Bartel DP. MicroRNA targeting specificity in mammals: determinants beyond seed pairing. Mol Cell. 2007;27(1):91–105.
Garcia DM, Baek D, Shin C, Bell GW, Grimson A, Bartel DP. Weak seed-pairing stability and high target-site abundance decrease the proficiency of lsy-6 and other microRNAs. Nat Struct Mol Biol. 2011;18(10):1139–46.
Morin RD, O’Connor MD, Griffith M, Kuchenbauer F, Delaney A, Prabhu AL, et al. Application of massively parallel sequencing to microRNA profiling and discovery in human embryonic stem cells. Genome Res. 2008;18(4):610–21.
Neilsen CT, Goodall GJ, Bracken CP. IsomiRs—the overlooked repertoire in the dynamic microRNAome. Trends Genet. 2012;28(11):544–9.
Heo I, Joo C, Cho J, Ha M, Han J, Kim VN. Lin28 mediates the terminal uridylation of let-7 precursor MicroRNA. Mol Cell. 2008;32(2):276–84.
Chang HM, Triboulet R, Thornton JE, Gregory RI. A role for the Perlman syndrome exonuclease Dis3l2 in the Lin28-let-7 pathway. Nature. 2013;497(7448):244–8.
Katoh T, Hojo H, Suzuki T: Destabilization of microRNAs in human cells by 3’ deadenylation mediated by PARN and CUGBP1. Nucleic Acids Res 2015.
Katoh T, Sakaguchi Y, Miyauchi K, Suzuki T, Kashiwabara S, Baba T, et al. Selective stabilization of mammalian microRNAs by 3’ adenylation mediated by the cytoplasmic poly(A) polymerase GLD-2. Genes Dev. 2009;23(4):433–8.
Jones MR, Blahna MT, Kozlowski E, Matsuura KY, Ferrari JD, Morris SA, et al. Zcchc11 uridylates mature miRNAs to enhance neonatal IGF-1 expression, growth, and survival. PLoS Genet. 2012;8(11), e1003105.
Burroughs AM, Ando Y, de Hoon MJ, Tomaru Y, Nishibu T, Ukekawa R, et al. A comprehensive survey of 3’ animal miRNA modification events and a possible role for 3’ adenylation in modulating miRNA targeting effectiveness. Genome Res. 2010;20(10):1398–410.
Koppers-Lalic D, Hackenberg M, Bijnsdorp IV, van Eijndhoven MA, Sadek P, Sie D, et al. Nontemplated nucleotide additions distinguish the small RNA composition in cells from exosomes. Cell Rep. 2014;8(6):1649–58.
Thornton JE, Du P, Jing L, Sjekloca L, Lin S, Grossi E, et al. Selective microRNA uridylation by Zcchc6 (TUT7) and Zcchc11 (TUT4). Nucleic Acids Res. 2014;42(18):11777–91.
Ha M, Kim VN. Regulation of microRNA biogenesis. Nat Rev Mol Cell Biol. 2014;15(8):509–24.
Hutvagner G. Small RNA asymmetry in RNAi: function in RISC assembly and gene regulation. FEBS Lett. 2005;579(26):5850–7.
Khvorova A, Reynolds A, Jayasena SD. Functional siRNAs and miRNAs exhibit strand bias. Cell. 2003;115(2):209–16.
Schwarz DS, Hutvágner G, Du T, Xu Z, Aronin N, Zamore PD. Asymmetry in the assembly of the RNAi enzyme complex. Cell. 2003;115(2):199–208.
Suzuki HI, Katsura A, Yasuda T, Ueno T, Mano H, Sugimoto K, et al. Small-RNA asymmetry is directly driven by mammalian Argonautes. Nat Struct Mol Biol. 2015;22(7):512–21.
Ahmed F, Ansari HR, Raghava GP. Prediction of guide strand of microRNAs from its sequence and secondary structure. BMC Bioinf. 2009;10:105.
Hu HY, Yan Z, Xu Y, Hu H, Menzel C, Zhou YH, et al. Sequence features associated with microRNA strand selection in humans and flies. BMC Genomics. 2009;10:413.
Frank F, Sonenberg N, Nagar B. Structural basis for 5’-nucleotide base-specific recognition of guide RNA by human AGO2. Nature. 2010;465(7299):818–22.
Huang DP, Ho JH, Poon YF, Chew EC, Saw D, Lui M, et al. Establishment of a cell line (NPC/HK1) from a differentiated squamous carcinoma of the nasopharynx. Int J Cancer. 1980;26(2):127–32.
Ressing ME, Keating SE, van Leeuwen D, Koppers-Lalic D, Pappworth IY, Wiertz EJ, et al. Impaired transporter associated with antigen processing-dependent peptide transport during productive EBV infection. J Immunol. 2005;174(11):6829–38.
Wildeman MA, Novalic Z, Verkuijlen SA, Juwana H, Huitema AD, Tan IB, et al. Cytolytic virus activation therapy for Epstein-Barr virus-driven tumors. Clin Cancer Res. 2012;18(18):5061–70.
Barturen G, Rueda A, Hamberg M, Alganza A, Lebron R, Kotsyfakis M, Shi B-J, Koppers-Lalic D, Hackenberg M: sRNAbench: profiling of small RNAs and its sequence variants in single or multi-species high-throughput experiments. In: Methods in Next Generation Sequencing. 2014;1:21-31.
Skalsky RL, Kang D, Linnstaedt SD, Cullen BR. Evolutionary conservation of primate lymphocryptovirus microRNA targets. J Virol. 2014;88(3):1617–35.
Gottwein E, Corcoran DL, Mukherjee N, Skalsky RL, Hafner M, Nusbaum JD, et al. Viral microRNA targetome of KSHV-infected primary effusion lymphoma cell lines. Cell Host Microbe. 2011;10(5):515–26.
We thank Lione Willems and Thierry Noorbergen for assistance in molecular cloning.
This work was supported by Veni grant 916.10.138 from The Netherlands Organization for Scientific Research to RJL and Dutch Cancer Society grant UU 2012-5667 to EJHJW and RJL. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Availability of data and materials
The data sets supporting the results of this article are available in the NCBI Sequence Read Archive repository, http://www.ncbi.nlm.nih.gov/sra, SRA study accession SRP066099, Bioproject PRJNA301879.
MH, RJL and EW conceived of the study, and participated in its design and coordination. EK and MH performed the experiments. MH and EK analysed the flow cytometry data. MH performed the analysis of the deep sequencing data. MH, RJL and EW wrote the paper. All authors read and approved the final manuscript.
The authors declare that they have no competing interests.
Consent for publication
Ethics approval and consent to participate
List of miRNA sequences and perfect target sequences. List of miRNA sequences annotated in miRBase or (*) based on deep sequencing data from Chen et al. 2010 (column 2). Perfect target sequences (column 3) were cloned downstream of the mCherry reporter gene to monitor miRNA expression. # The deep-sequencing analysis in Fig. 3 identified a single nucleotide change (position 17, T changed to C) in the mature miR-BART19-5p miRNA sequence in the C666-1 strain. Therefore, the miR-BART19-5p miRNA sensor used in the C666-1 cells to monitor miRNA activity does not fully match the miR-BART19-5p sequence expressed in these cells (DOCX 20 kb)
Good reproducibility between miRNA reporter assays. Shown are the log-transformed relative reporter expression values of two independent experiments (each representing three technical replicates) in each of the four EBV+ cell lines. The data points are located close to the diagonal dashed lines, indicating good reproducibility (PDF 121 kb)
Comparison of functional EBV miRNA expression between AKBM and other EBV+ cell lines. miRNA reporter expression (data also shown in Fig. 1) was compared between AKBM and the three other EBV+ cell lines. Every dot represents the relative reporter expression of a specific miRNA as a percentage of the control reporter. The dashed line indicates the diagonal (no difference in relative reporter expression between both cell lines) (PDF 91 kb)
Assessment of miRNA reporters in EBV− cells. Reporter activity for miR-BART2-3p, miR-BART3-5p and controls was compared between AKBM and the EBV− Burkitt’s lymphoma cell lines, 2A8 and AK31 (AKBM data also shown in Fig. 1). The miR-BART2-3p reporter is only downregulated in AKBM but not in the other two cell lines, suggesting that downregulation of the reporter is specific for EBV+ cells and thus likely due to miR-BART2-3p expression. The miR-BART3-5p reporter was downregulated in all three cell lines, suggesting the influence of a common host cell factor such as a human miRNA on reporter expression (PDF 101 kb)
List of EBV miRNA sequence variants. The most highly abundant sequence, the normalized sequence read counts (reads per million genome-mapped reads), the percentage of canonical sequences and percentage of 5’ terminus variants of each EBV miRNA in the five cell lines are shown (XLSX 22 kb)
AKBM cells are EBV positive. A) Schematic representation of the EBV super-sensor. Five perfect EBV miRNA target sites were cloned downstream of the mCherry reporter gene. The vector allows for sensing of EBV presence by monitoring the combined functional expression of 5 EBV miRNAs. B) EBV-positive AKBM cells and EBV negative AK31 and 2A8 cells were either left untreated (cells alone), were transduced with an mCherry-control vector (control sensor), or were transduced with the EBV sensor from a). Cells were monitored for mCherry reporter expression at 7 days post transduction by flow cytometry (PDF 165 kb)
Comparison of deep sequencing reads between different cell lines. Deep sequencing reads were compared between the four EBV+ cell lines and the EBV− cell line HK1-miR-BART. Every dot represents a single miRNA arm. BHRF1 miRNAs, BART2-3p and BART2-5p are not depicted in the HK1-miR-BART plots as these are not expressed in these cells. Pearson correlations R square and p-values are presented. Good correlations were observed between the read counts for BART miRNAs in all combinations of cell lines (PDF 157 kb)
About this article
Cite this article
Hooykaas, M.J.G., Kruse, E., Wiertz, E.J.H.J. et al. Comprehensive profiling of functional Epstein-Barr virus miRNA expression in human cell lines. BMC Genomics 17, 644 (2016). https://doi.org/10.1186/s12864-016-2978-6
- Epstein-Barr virus (EBV)
- miRNA reporter
- miRNA sensor
- Small RNA sequencing