1 Introduction to consensusDE

consensusDE aims to make first pass differential expression (DE) analysis, with reporting of significance scores from multiple methods easy. It implements wrappers for Voom, DEseq2 and EdgeR and reports differential expression results seperately, as well as merging the results into a single table for determining consensus. The results of the merged table, are ordered by the summed ranks of the p-values for each algorithm.

Core functionality is simplified into two function:

  • buildSummarized() generate a summarized experiment that counts reads mapped (from bam files) against a transcriptome
  • multi_de_pairs() perform DE analysis (all possible pairwise comparisons)

Below, we demonstrate the core functionality of consensusDE as well as how to plot results from obtained results using the diag_plots() function.

2 consensusDE examples

Begin by first installing and then loading the consensusDE library. To illustrate functionality of consensusDE, we will utilise data from the airway and annotation libraries as follows. Begin by installing and attaching data from these libraries:

3 Building a summarized experiment

A summarized experiment is an object format that stores all the relevant information for performing differential expression analysis. buildSummarized() allows users to build a summarized object by simply providing 1) a table of bam files (more below on format), 2) a directory of where to locate the bam files and 3) a transcript database to map the reads to (either a gtf file or txdb). We will use bam files attached to this package (from GenomicAlignments) as an example:

The minimum information is now ready to build a summarized experiment as follows:

This will output a summarized object that has mapped the reads for the bam files that are listed in sample_table, located in bam_dir, against the transcript database provided: TxDb.Dmelanogaster.UCSC.dm3.ensGene. Bam file format, whether “paired” or “single” end (the type of sequencing technology used) must be specified using the read_format parameter. gtf formatted transcript databases can also be used instead of a txdb, by providing the full path to the gtf file using the gtf parameter. To save the summarized experiment externally, for future use, specify the path to save the summarized experiment using output_log To see details of all parameters see ?buildSummarized.

Overview of the summarized experiment:

3.1 Filtering low count data

buildSummarized() also allows users to filter out low read counts. This can be done when building the summarized experiment, or re-running with the summarized experiment output using buildSummarized(). See “Performing Differential Expresssion” below with filter example.

4 Performing Differential Expresssion

For differential expression (DE) analysis we will use the airway data for demonstration. See ?airway for more details for this experiment. NOTE: the summarized meta-data must include the columns “group” and “file” to build the correct models. For illustration, we sample 1000 genes from this dataset.

Running multi_de_pairs() will perform DE analysis on all possible pairs of “groups” and save these results as a simple list of “merged” - being the merged results of “deseq”, “voom” and “edger”, as well as the latter three as objects independently. To access the merged results:

4.1 Annotating DE results

It is often useful to add additional annotated information to the output tables. This can be acheived by providing a database for annotations via ``ensembl_annotate. Annotations needs to be a Genome Wide Annotation object, e.g. org.Mm.eg.db for mouse or org.Hs.eg.db for human from BioConductor. For example, to install the database for the mouse annotation, go to http://bioconductor.org/packages/org.Mm.eg.db and follow the instructions. Ensure that after installing the database package that the library is loaded usinglibrary(org.Mm.eg.db). When running, "'select()' returned 1:many mapping between keys and columns" will appear on the command line. This is the result of multiple mapped transcript ID to Annotations. Only the first annotation is reported. See?multi_de_pairs``` for additional documentation.

An example of annotating the above filtered airway data is provided below:

4.2 Writing tables to an output directory

multi_de_pairs provides options to automatically write all results to output directories when a full path is provided. Which results are output depends on which directories are provided. Full paths provided to the parameters of output_voom, output_edger, output_deseq and output_combined will output Voom, EdgeR, DEseq and the merged results to the directories provided, respectively.

5 Removing unwanted variation (RUV)

consensusDE also provides the option to remove batch effects through RUVseq functionality. consensusDE currently implements RUVr which models a first pass generalised linear model (GLM) using EdgeR and obtaining residuals for incorporation into the SummarizedExperiment object for inclusion in the models for DE analysis. The following example, uses RUV to identify these residuals. To view the residuals in the model see the resisuals section below in the plotting functions. Note, that if ruv_correct = TRUE and a path to a plot_dir is provided, diagnostic plots before and after RUV correction will be produced. The residuals can also be accessed in the summarizedExperiment as below. These are present in the “W_1” column. At present only one factor of variation is determined.

6 Plotting functions

When performing DE analysis, a series of plots (currently 10) can be generated and saved as .pdf files in a plot directory provided to multi_de_pairs() with the parameter: plot_dir = "/path/to/save/pdfs/. See ?multi_de_pairs for description.

In addition, each of the 10 plots can be plotted individually using the diag_plots function. See ?diag_plots for description, which provides wrappers for 10 different plots. Next we will plot each of these using the example data.

6.1 Mapped reads

Plot the number of reads that mapped to the transcriptome of each sample. The sample numbers on the x-axis correspond to the sample row number in the summarizedExperiment built, accessible using colData(airway). Samples are coloured by their “group”.

6.4 RUV residuals

Residuals for the RUV model can be plotted as follows: