MeT-DB V2.0: elucidating context-specific functions of N6-methyl-adenosine methyltranscriptome.

Methyltranscriptome is an exciting new area that studies the mechanisms and functions of methylation in transcripts. A knowledge base with the systematic collection and curation of context specific transcriptome-wide methylations is critical for elucidating their biological functions as well as for developing bioinformatics tools. Since its inception in 2014, the Met-DB (Liu, H., Flores, M.A., Meng, J., Zhang, L., Zhao, X., Rao, M.K., Chen, Y. and Huang, Y. (2015) MeT-DB: a database of transcriptome methylation in mammalian cells. Nucleic Acids Res., 43, D197-D203), has become an important resource for methyltranscriptome, especially in the N6-methyl-adenosine (m6A) research community. Here, we report Met-DB v2.0, the significantly improved second version of Met-DB, which is entirely redesigned to focus more on elucidating context-specific m6A functions. Met-DB v2.0 has a major increase in context-specific m6A peaks and single-base sites predicted from 185 samples for 7 species from 26 independent studies. Moreover, it is also integrated with a new database for targets of m6A readers, erasers and writers and expanded with more collections of functional data. The redesigned Met-DB v2.0 web interface and genome browser provide more friendly, powerful, and informative ways to query and visualize the data. More importantly, MeT-DB v2.0 offers for the first time a series of tools specifically designed for understanding m6A functions. Met-DB V2.0 will be a valuable resource for m6A methyltranscriptome research. The Met-DB V2.0 database is available at http://compgenomics.utsa.edu/MeTDB/ and http://www.xjtlu.edu.cn/metdb2.

INTRODUCTION

Methyltranscriptome is an exciting, emerging area that studies methylation in the transcriptome. In contrast to well-established DNA methylation, transcriptome methylation is largely an uncharted territory. Among different types, N6-methyl-adenosine (m6A) is the most abundant and intensively studied transcriptome methylation, existing in transcriptomes of mammalian and other organisms. Through methylated RNA immunoprecipitation sequencing (MeRIP-seq) (1), m6A has been found to exist in >25% mRNAs in mammalian cells, particularly enriched near stop codon and with a consensus RRACH motif (R = G or A; H = A, C or U) (2). m6A is also found to be highly dynamic (3). It is catalyzed by ‘writers’ or m6A methylases METTL3 and METTL14, two interacting with the Wilms Tumor 1 Associated Protein (WTAP) (4) to form methyltransferase complex (5). This methylation can also be reversed by ‘erasers’ including the alkylated DNA repair protein (ALKBH5) (6) and the obesity associated protein (FTO) (3) (less effective in vivo (7)). In addition to writers and erasers, m6A-binding proteins such as YTH protein family have been identified as m6A readers, which can mediate the biological function of m6A through selectively recognizing and binding to m6A. A close relationship between m6A and mRNA metabolism has also been established, and m6A is involved in regulating diseases and virus infection (5,6,8–14). These exciting findings have spurred intense research recently on roles of m6A in regulating different physiological processes as well as their potential as therapeutic targets (5,6,11,14–21).

MeT-DB (22), established in 2014, was the first comprehensive database focusing on m6A methyltranscriptome. Since then, a large additional number of MeRIP-seq datasets produced under different experimental conditions have been released (Supplementary Figure S1). Moreover, a series of bioinformatics tools have been developed for predicting m6A peaks (exomePeak (23), MeTPeak (24)) and differential m6A analysis (exomePeak, MeTDiff (25)), for visualizing the characteristics of peaks in transcripts (Guitar (26)), and for predicting context-specific m6A driver genes and networks (m6A-Driver (27)).

Here, we present the second version of MeT-DB, which has been entirely redesigned to focus more on helping elucidate context-specific m6A functions (Figure 1). Compared with MeT-DB, MeT-DB V2.0 has ∼2.5× MeRIP-seq samples and ∼7.6× predicted m6A peaks from 7 species and 26 studies. Particularly, a new database, TREW, which includes over 118k targets of eight different m6A readers, erasers and writers is also integrated. MeT-DB V2.0 also expands the collections of other functional data such as micro-RNA target sites, Single nucleotide polymorphisms (SNPs), binding sites of splicing factor as well as RNA-binding proteins (RBPs), and information about cancer genes. MeT-DB v2.0 adopts a table view interface with multiple query options to deliver diverse information about m6A in parallel. MeT-DB v2.0 also replaces the original genome browser with a more efficient and powerful genome browser to be able to display 979 tracks for all species in a standard manner. More importantly, MeT-DB v2.0 also offers for the first time a series of tools specifically designed for understanding m6A functions. We discuss next the detailed improvements in database and web interface.

Figure 1.

Overall design of MeT-DB V2.0 database. MeT-DB V2.0 is composed of the database and web interface. The MeT-DB database includes the core database that contains context-specific m6A peaks and single-base sites, the TREW database that contains target sites of m6A readers, writers and erasers, and the functional database such as micro-RNA target sites, binding sites of RNA binding proteins and information about cancer genes. There are three functional modules in the web interface: table view facilitates researcher to explore and search the data in detail, the genome browser helps the user visualize and compare m6A peaks and functions data, and the tool module includes two useful web servers for investigating the functions of m6A methyltranscriptome.

MATERIALS AND METHODS

MeRIP-seq data processing

MeRIP-seq experimental information of all collected studies was obtained from the original papers or NCBI Gene Expression Omnibus, while raw sequencing data samples were downloaded from short read archive.

To detect m6A peaks, sequencing data quality was first evaluated by FASTQC (v0.11.4). Adaptors or low quality nucleotides were removed by Trim Galore (v0.4.2) according to the evaluation results of FastQC. Then, reads in the IP/Input FASTQ files were aligned to the genome by Tophat2 (v2.1.0) (28) with default options to generate IP/Input BAM files. BAM files were subsequently converted to bigwig files for visualization. Peak calling was performed on the input and IP BAM files by exomePeak (23). For each predicted m6A peak, its chromosomal location including start/end position, strand information, P-value, fold enrichment and q-value (FDR) were reported. For each sample, sequence motifs of predicted m6A peaks were obtained using the MEME (v4.11.2) (29) suite and the peak distribution at a transcript level was also plotted by the Guitar package (26).

Single-base m6A sites were also predicted by searching the RRACH motif in peaks identified by exomePeak. Transcript sequences of the peak region that contain only exons were first extracted, from which the location of RRACH motifs was identified. The genome positions of ‘A’ in the identified motifs were annotated as single-base m6A sites. The distances of the predicted single-base m6A sites to their corresponding peak center were also calculated as the confidence scores of the prediction.

Transcriptome-wide expression levels for each sample were also calculated based on the aligned BAM files of the MeRIP-seq input samples; cufflinks (v4.11.2) (30) with default settings was employed to calculate the gene/isoform expression Fragments Per Kilobase of transcript per Million mapped reads (FPKM) values and the reads counts generated by HTSeq (v0.6.1) (31) were also provided to facilitate further analysis.

TREW data preparing

TREW or the Target of m6A readers, erasers and writers is our newly constructed database about the binding sites of m6A methyltransferases (METTL3, WTAP, METTL14 and KIAA1429), demethylases (FTO and ALKBH5) and readers (YTH family proteins). (Supplementary Figure S2). To determine the target sits, ParCLIP-seq data were retrieved directly from original publications, where the raw data were first processed with Trim Galore and FASTX-Toolkit (v0.0.13) for quality control, and then aligned to human hg19 or mouse mm10 reference genome respectively with Tophat2. Also, differential m6A analysis was performed with exomePeak and QNB (32) packages under the default setting on MeRIP-seq data of m6A methylase or demethylase perturbation. The significant differential m6A peaks after perturbation were determined to the target peaks.

Techniques

All datasets were processed and stored in a MySQL Database Management System installed on an X86–64 server with CentOS Linux OS. The database consists of 60 tables that comprise ∼41 million records. The database query, genome browser and user interface were developed using PHP, JavaScript, jQuery and Bootstrap frameworks.

RESULTS

Content of database

Core MeT-DB V2.0 database

We expanded the number of MeRIP-seq samples in this update from 74 to 185, which come from 26 independent studies covering 7 species (Supplementary Figure S3 and Table S1) (1,2,9,11,13–15,17,19–21,33–47). More than 2.6 million m6A peaks and 1.1 million single-base sites were identified from these 185 samples. Annotation information was generated, including the transcript identifiers (UCSC ID, Entrez Gene ID, Gene Symbol and RefSeq ID), the transcript location information (5′UTR/CDS/3′UTR) and information regarding if the peak is on mRNA/lncRNA. Despite considerable variation in the numbers of detected m6A peaks across different MeRIP-seq conditions, the most enriched consensus sequences are all similar to the previously reported RRACH motif of m6A, suggesting that the exomePeak detected m6A peaks should have high specificity.

All the processing results, including m6A peaks, single-base m6A sites, motif, peak distribution plot, gene expression profiles, were stored in the MeT-DB core database and can be downloaded from the web page. The scale of data with comparison with MeT-DB V1.0 is summarized in Table 1.

Table 1.

Total number of identified m6A peaks and single-base m6A sites

	Included		Studies		Samples		m6A Peaks		Single-base m6A in V2.0
Species	V1.0	V2.0	V1.0	V2.0	V1.0	V2.0	V1.0	V2.0	Total	AVG. supp1	Genome assembly
Homo sapiens	Y	Y	4	13	16	75	2.1E + 5	1.3E + 6	4.3E + 5	14.5	hg19/hg38
Mus musculus	Y	Y	2	8	6	52	1.4E + 5	1.1E + 6	5.4E + 5	10.7	mm9/mm10
Sus scrofa		Y		1		3		1.3E + 4	2.9E + 4	2.4	susScr3
Danio rerio		Y		1		5		2.4E + 4	3.7E + 4	2.8	danRer10
Drosophila melanogaster		Y		1		2		83	92	1.4	dm6
Saccharomyces cerevisiae		Y		1		34		1.8E + 5	6.9E + 4	11.4	R64–1-1
Arabidopsis thaliana		Y		3		14		7.9E + 4	7.4E + 4	3.5	tair10
Total	2	7	5	26	22	185	3.5E + 5	2.6E + 6	1.1E + 6

Comparison of MeRIP-seq data collection and procession results between MeT-DB V1.0 and MeT-DB V2.0.

1AVG. supp (average support number) stands for average number of MeRIP-Seq samples that support each single-based m6A site.

TREW database

We collected ParCLIP-seq and MeRIP-seq samples for eight regulator/reader proteins (including FTO, KIAA1429, METTL14, METTL3, WTAP, HNRNPC, YTHDC1 and YTHDF1) from 10 independent studies (Supplementary Table S2) (8,9,15,20,21,36,47–50). Finally, a total of 118, 164 m6A modifications sites that preferentially targeted by a specific protein regulator are collected or detected. The target sites are further annotated with transcript regions (5′UTR, CDS, 3′UTR, stop codon, transcription start sites and miRNA target site) and RNA types (tRNA, mRNA, lncRNA and sncRNA). Target sites were also converted to hg38 and mm9 by UCSC liftOver tool (51). This information is essential for understanding the biological functions of m6A.

Functional data

To help understand the regulatory roles of m6A, we also integrated the following six relevant functional datasets into the Met-DB v2.0 database.

miRNA target sites

Predicted miRNA target sites from TargetScan (version 7.1) (52) for human and mouse as well as from miRanda for human, mouse and fly (August 2010 release) (53) were included. Furthermore, experimentally validated miRNA and target genes interaction pair information for human, mouse, zebrafish and fly were downloaded from miRTarBase (54).

Splicing factor binding sites

A total of 655 and 125 binding sites of human and mouse splicing factors, respectively, were obtained from SpliceAidF(v1.1 03/2013) (55). Each site includes the name of the binding splicing factor and the genome location of the binding site.

SNP. 40627 human literature-derived collected SNP-trait associations of 30044 SNPs obtained from GWAS (All associations v1.0.1) (56) were included.

RBP-binding site

PAR-CLIP and HITS-CLIP detected binding locations of 24 human and 5 mouse RBPs were retrieved from StarBase version 2.0 (57) and included in the database.

Cancer related genes

A total of 761 human and 628 mouse tumor suppressor genes were obtained from TSGene database version 2.0 (58). Also, 576 cancer genes were downloaded from COSMIC (v82) (59).

These relevant functional data are displayed as part of the query output of the core Met-DB database. Each functional dataset can also be downloaded from the website.

Table view interface and query the database

To facilitate the examination of the context-specific samples and the query of the entries of the database, MeT-DB V2.0 included a table view interface (Supplementary Figure S4) based on DataTables jQuery plug-in. A table view presents information of every genomic features of an m6A peak/site as a row, thus simultaneously offering a large amount of information about m6A. It can also show the potential interactions between a particular m6A peak/site with other transcriptome features and functional data. Presenting this rich information can provide multifaceted perspectives of m6A methylation and help generate hypotheses of their potential biological functions. Moreover, the table view provides flexible ways to search the large database, making Met-DB V.20 a powerful platform for discovering the biological functions of m6A.

Three different ways of querying the database are made available through the MeT-DB web page, namely by samples, by m6A peaks and by single-base m6A sites, providing information about m6A methyltranscriptome at different scales ranging from a global perspective to a single-nucleotide resolution. The search-by-sample function aims to provide the user with a comprehensive view of context-specific m6A methyltranscriptome. The query takes sample ID, experiments, species, cell line or tissue as input and returns transcriptome-wide information about m6A. For each returned sample, the user can investigate the m6A peak distribution at a transcript level and m6A peak sequence motif and download the BED file, detailed information of the transcriptome-wide predicted m6A peaks, the gene and isoform expression FPKMs and sample reads counts. The search-by-peak and search-by-site functions instead deliver information about predicted m6A peaks and single-base site integrated from all the samples and display it in the table view interface. For each peak/site, the genomic location, associated gene name/ID, peak enrichment fold change, prediction confidence and the number of other peaks/sites from other samples that overlap with this queried m6A peak/site are displayed in the columns of the table. To further facilitate functional discovery, additional functional information including overlapped binding sites of m6A readers, writer and erasers, binding sties of RNA-binding proteins and slicing factors, miRNA target sites, SNPs and status of its association with a tumor suppressor or a cancer gene is also provided in the columns of the table view. Column specific search, located at the bottom of the table, further enables the search within a corresponding column. Advanced filtering functions are also implemented to allow the user to further restrict their search by any entry attribute such as genome location, P-value, score, gene name, transcript name, etc. These query functions should allow users to screen out most relevant elements from the huge information stored in the MeT-DB database. Users can also export data by clicking export buttons under any search conditions for offline investigation. Besides, users can view the detailed information of a specific entry in the genome browser by clicking the genome browser icon at the very left end of each row.

Genome browser

Met-DB originally included a self-constructed genome browser to graphically illustrate the genomic features and context of m6A peaks. To provide an efficient service on MeT-DB data of a much larger scale, Met-DB v2.0 replaced the old genome browser with JBrowse (60). JBrowse is a lightweight open source, JavaScript-based genome browser that effectively renders most genomic features on the client machine rather than the server. Besides, JBrowse contains more dynamic visualization capabilities and commonly used functions. With this significant update, the response speed of genome browser is much quick comparing to the old version.

All information of MeRIP-Seq data, TREW and other functional datasets have been converted to BED, GFF or BigWig format and imported into JBrowse as tracks. Furthermore, categories as well as attributes of each track have been described explicitly for generating track filter function module. Based on this module, users can screen out their interested ones among hundreds of tracks by JBrowse track filtering panel.

Guitar plot

The ability to visualize the distribution of m6A at a transcript level and compare the differences of such distribution under the different condition is very important in generating new hypotheses about m6A functions. Guitar plot (26) is a tool designed for visualizing the transcript level m6A distribution. It can also be used to visualize distributions of any other type of genome features and transcriptome methylations stored in BED files. MeT-DB V2.0 includes a web server for Guitar to generate the plots of m6A distribution from user custom data. The web server brings the full capability of Guitar to generate not only the distribution plot of a single sample but also the plot that compares distributions from multiple samples. The plot is in PDF format and is publication ready.

m6A-Driver

m6A-Driver (27) is a published tool for predicting m6A-driven genes and associated networks, whose functional interactions are likely to be actively modulated by m6A methylation under a particular condition. MeT-DB V2.0 also includes a web server for m6A-Driver to generate m6A-driven genes and associated networks from user custom data samples. The inputs to the m6A-Driver web server are txt files, each containing a list of official gene symbols, representing a replicate sample of context-specific m6A targeted genes under a case-control condition obtained from, for instance, differential m6A analysis using exomePeak or MetDiff. m6A-Driver predicts m6A driven genes from the provided gene lists by assessing their topological and biological significance using a random walk with restart algorithm applied to the protein–protein interaction network. The output of m6A-Driver web server contains two files; one is a text that contains m6A driven genes and the edges of underlying m6A driven gene network and the second file is the figure of the network.

Data download and export

All MeT-DB V2.0 data samples including transcriptome-wide m6A peaks and single-base sites, TREW dataset and six functional data files are available for downloading on the download web page. Transcript-level peak distribution, m6A sequence motif, gene/isoform FPKM expression levels, reads counts for each MeRIP-seq experiment can be downloaded in detail information region of MeT-DB Core experiments table view web page. User query or filtering results of each table can also be downloaded in a comma-separated value or tab-delimited text format, accessible via the export buttons above the corresponding table.

CONCLUSION

MeT-DB V2.0 is a comprehensive and significantly enhanced database with a redesigned web interface. By collecting and integrating more MeRIP-seq samples, functional datasets and tools in one place, we believe that MeT-DB V2.0 can help accelerate the discovery of the unrecognized regulatory roles of m6A and foster the further development of bioinformatics tools for methyltranscriptome research. Comparing with the RMBase (61), another database for RNA modifications, MeT-DB V2.0 not only included more MeRIP-seq experiments from different studies, species and conditions, but also provided more functional data, additional web tools and convenient access to condition-specific RNA methylation information for functional studies.

m6A and methyltranscriptome research is a fast evolving area, with exciting discoveries and new datasets being published with an unprecedented pace. To appropriately reflect the development of this research, the MeT-DB database needs to be continuously updated and improved. Specific areas that will see further updates include: (i) continuously collecting published MeRIP-seq samples and carrying out predictions of m6A peaks/sites, (ii) including other types of post-transcriptional methylations once the corresponding data samples are enough for us to build a high-quality database, (iii) improving user interface to improve the ability to mine information from database and (iv) providing more web server tools for methyltrancriptome data analysis and functional prediction. We envision that MeT-DB would become the central resource for data and discovery hub for methyltrancriptome research.

AVAILABILITY

The Met-DB V2.0 database is available at http://compgenomics.utsa.edu/MeTDB/ and http://www.xjtlu.edu.cn/metdb2.

Click here for additional data file.