RNA editing is the molecular process by which ribonucleic acid (RNA) molecule’s nucleotide sequences are modified enzymatically after RNA has been generated by RNA polymerase. There are two types of canonical RNA editing: adenine-to-inosine (A-to-I) and cytosine-to-uracil (C-to-U). The mechanisms of these two types of RNA editing have been thoroughly studied. The A-to-I editing is mediated by adenosine deaminase acting on RNA (ADAR) (Bass and Weintraub, 1987), and C-to-U editing happens when cytidine deaminase deaminates a cytidine base into a uridine base (Gott and Emeson, 2000). A-to-I RNA editing accounts for >95% of all reported RNA editing events and over four hundred million adenosine sites in the human genome are presumed to be subject to A-to-I editing (Bazak et al., 2014; Picardi et al., 2017). The mapping of whole-genome RNA editing events leads to the elucidation of inosinomes (Picardi et al., 2015) or editomes (Diroma et al., 2017).

Most well characterized functional editing sites are found in transcripts for neuronal transporters and channel proteins in the brain (Meier et al., 2016). Recent advances in RNA editing research have revealed functional effects of A-to-I RNA editing in extended contexts, especially tumors. Peng et al. (2018), demonstrated experimentally that non-synonymous A-to-I RNA editing can result in alternative protein sequences by modifying amino acids in cancer, and may subsequently affect drug sensitivity regulation (Han et al., 2015). Furthermore, increased editing activity has been associated with poor cancer prognosis (Paz-Yaacov et al., 2015). A comprehensive review of the roles of RNA editing in human cancers was published recently (Qian et al., 2018).

Given a large amount of paired DNA-RNA sequencing data accumulated over the last decade, it became feasible to detect RNA-DNA Difference from direct comparison of DNA and RNA nucleotide sequences, and millions of RDDs have been documented in RNA editing databases, such as REDIportal (Picardi et al., 2017). REDIportal was based on the Genotype-Tissue Expression (GTEx) data (Lonsdale et al., 2013) which has been used for comprehensive RNA editing analysis (Tan et al., 2017). Other than studying the occurrence of RNA editing events (qualitative identification), there is an increasing interest in examining the level of RNA editing, especially in experiments that involve multiple conditions (quantitative comparison). Technically speaking, the “RNA editing level” of an RNA sample can be calculated as the percentage of edited bases at an RNA editing site, which is an index specific to a particular biological condition. For instance, JACUSA (Piechotta et al., 2017) was empowered to compare RNA editing level across different experimental conditions, thus detecting RNA-RNA differences. In a study of Wood-decaying fungi, JACUSA was employed to compare RNA editing levels across varied culture conditions and it helped confirm the more prominent impact of culture conditions than host wood species (Peng et al., 2018). In analogy to the pervasive differential expression analyses in transcriptome studies, the rising editome studies will soon be pivoting on differential editing analyses.

Abundant high-throughput sequencing data of human cancers enable us to study the differential RNA editing level in tumor samples relative to normal. Here, in various perspectives, we performed RNA editing level analyses on 6,236 samples from 14 cancers from TCGA, aiming to obtain global views on RNA editing levels across different sample sources (tumor vs. health), genomic locations (exonic, intronic, and intergenic), and gene features (protein-coding genes, lincRNAs, and pseudogenes). Because of the dominance of A-to-I RNA editing, our study primarily focused on the canonical A-to-I editing type.

We obtained 230,079 RNA editing events across 107,049 genomic positions from a previous RNA editing study by Han et al. (2015). These RNA editing sites were identified from 6,236 samples of 14 cancers from TCGA. These samples were collected from 5,672 subjects, 544 of which had paired tumor and normal samples. The RNA editing detection methodology and filters were described in detail by Han et al. (2015). The RNA editing sites obtained from Han et al. (2015) study were the sites passed their quality control measures, and used to generate their results. Our study is an extension of Han et al. (2015) study, while verifying some of the previous conclusions, we focused more on the non-coding regions.

The level and rate of RNA editing were computed using the following formulas:

Gene edit frequency was defined as number of RNA edit within the gene divided by the gene length. Hyper RNA edited genes were identified by rank the gene edit frequency from large to small.

Furthermore, 4.66 million known RNA editing events from database REDiportal (Picardi et al., 2017) were downloaded and summarized for comparison purpose. Built primarily upon a different source “GTEx project.” The detailed RNA editing identification method is described in REDiportal publication (Picardi et al., 2017).

In the present study, RNA editing in chromosome 1–22, X and Y were studied thoroughly, while the mitochondrial genome was omitted. Annotation of the RNA editing positions was done using ANNOVAR (Wang et al., 2010). Wilcoxon Rank Sum test was used to compare the RNA editing level between tumor and normal. A Benjamini Hochberg adjusted p-value <0.05 was used as the significant threshold.

Our analyses found that the RNA editing rate is highly associated with gene density, especially protein-coding gene density. Chromosome 19 was found to harbor the highest RNA editing rate among all chromosomes (excluding mitochondria). The highest RNA editing rate can be explained by the highest gene density in chromosome 19. The correlations between RNA editing rate and protein-coding gene density is higher than the correlations between RNA editing rate and lincRNA density. This could result from the bias of the source data. The TCGA RNA editing events were detected by comparing exome sequencing and RNA-seq data. Exome sequencing targets protein coding regions primarily. RNA-seq, depending on library preparation method, can target protein coding or total RNA. However, when using exome sequencing data as a reference, the detection range for RNA editing events is limited to protein-coding regions. Furthermore, the RNA editing detection is also limited to gene expression. No RNA editing can be detected with non-expressed genes. Frankly, the majority of the RNA editing source sequencing data were not designed for RNA editing detection initially but were data mined for RNA editing at a later time. The focus of protein-coding regions of the source data resulted that the current detected RNA editing events heavily reside in protein coding regions.

A majority of the RNA editing events were found in either intronic regions or 3′-UTRs for TCGA and GTEx data. This ubiquitous phenomenon has been discovered in previous findings (Alseth et al., 2014; Daniel et al., 2015). Significant proportional difference was observed for RNA editing in 3′-UTR regions between TCGA and GTEx, also suggesting potential difference of RNA editing pattern between cancerous and healthy subjects. However, such difference can also arise from a potential limitation of GTEx data. GTEx derived genomic data from tissues acquired via autopsy of donors. The GTEx inclusion criteria required tissue collection within 24 h of death. Reports (Broniscer et al., 2010; Birdsill et al., 2011) have shown that post mortem interval is negatively correlated with RNA integrity and the degradation of RNA can lead to loss of the RNA expression signal (Durrenberger et al., 2010). The effect of post-mortem interval from GTEx on RNA editing has yet to be documented.

The intronic RNA editing is less likely to cause functional consequence. The 3′UTR contains key regulatory elements such as miRNA response elements, AU-rich elements, and poly(A) tails. These elements can affect mRNA stability, export, and translation efficiency. Thus, it has been proposed that RNA editing in 3′UTR can up-regulate mRNA expression by acting as a miRNA sponge, repressing miRNA expression (Borchert et al., 2009; Wang et al., 2013; Zhang et al., 2016). It has also been suggested that RNA editing in 3′UTR can stabilize mRNA by modulating the RNA secondary structure (Brummer et al., 2017). Our across tissue type analysis of RNA editing shows that there could be tissue exclusivity of RNA editing, ∼44% of the RNA editing positions occurred in one or two types of cancer. A small portion (3.4%) of the RNA editing positions are ubiquitous across all cancer types.

There have been multiple reports that RNA editing events tend to occur in close vicinity (Carmi et al., 2011; Danecek et al., 2012; Porath et al., 2014; John et al., 2017), these phenomena have been referred to as clusters, editing islands or hyper-editing regions. We found many of the protein-coding genes had high RNA editing rates after adjusting for the gene length, many of these RNA editing causes non-synonymous changes, which has been proven to generate alternative protein sequences (Peng et al., 2018) which can potentially significantly alter function, including drug sensitivity (Han et al., 2015). Notably, 9 out of the top 10 RNA editing protein coding genes from TCGA had prognostic value (Table 2). However, we cannot assert that the RNA editing modified gene expression. The exact mechanisms and potential functions of RNA editing in these protein coding genes remain to be discovered.

A portion of the RNA editing events do occur in non-coding RNAs, which have gained substantial interest in recent years. One of the common types of non-coding RNA that we explored in our study is lincRNA, which are traditionally thought to be non-functional. The increasing reports of lincRNA-disease associations over the last decades suggest that there are unknown mechanisms by which lincRNA exerts its influence (Wang et al., 2017). There have been suggestions that lincRNAs function through regulation of protein coding genes by trans expression quantitative trait loci (eQTL) (Cabili et al., 2011; Hon et al., 2017). Another type of non-coding RNA we examined was the pseudogenes. The known potential primary function of pseudogenes involve them acting as miRNA sponges, alleviating the regulation effect of miRNAs (Poliseno et al., 2010; Karreth et al., 2015). Our results show that lincRNAs and pseudogenes are among the top hyper RNA edited genes in GTEx. Whether RNA editing affects the potential functions of lincRNA and pseudogene has been understudied and will require further exploration.

This study uses publically available data downloaded from GDC data repository. It does not contain any patient information or patient data.

MC, YG, DS, and Y-yZ wrote the manuscript. MC, HY, and WY performed the analyses and designed the figures and tables for the first submission. YG designed the study. SN, QS, and OO performed the additional analyses and designed the new figures during revisions.

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.