Page Comparison

Overview

In this lab, we'll look at how to use cummeRbund to visualize our gene expression results from cuffdiff.  CummeRbund is part of the tuxedo pipeline and it is an R package that is capable of plotting the results generated by cuffdiff.

Image Added

What can you do with cummeRbund?

Get your data

cds
cd my_rnaseq_course
cp -r /corral-repl/utexas/BioITeam/rnaseq_course/cuffdiff_results .
 
 

#We need the cuffdiff results because that is the input to cummeRbund.

A) First load R and enter R environment

module#module load R will not work for us because we need the latest version of R. So, we'll load a local version.

$BI/bin/R

 
What is $BI here?

B) Within R environment, set up install package cummeRbund

source("http://bioconductor.org/biocLite.R")
biocLite("cummeRbund")

C) Load cummeRbund library and read in the differential expression results.  If you save and exit the R environment and return, these commands must be executed again.

library(cummeRbund)
cuff_data <- readCufflinks('diffcuffdiff_outresults')

D) Use cummeRbund to visualize the differential expression results.

...

And after executing the necessary commands, add the line:

dev.off()

...

Exercise 1: Generate a scatterplot comparing genes across conditions C1 vs C2.

pdf(file="scatterplot.pdf")

csScatter(genes(cuff_data), 'C1', 'C2')

dev.off()

Let's look at our scatterplot here.

Exercise 2:  Pull out from your data, significantly differentially expressed genes and isoforms.

gene_diff_data  <- diffData(genes(cuff_data))
sig_gene_data  <- subset(gene_diff_data, (significant ==  'yes'))

#Count how many we have
nrow(sig_gene_data)

isoform_diff_data <-diffData(isoforms(cuff_data))
sig_isoform_data <- subset(isoform_diff_data, (significant == 'yes'))

#Count how many we have
nrow(sig_isoform_data)

Exercise 3: For a gene, regucalcin, plot gene level expression.

pdf(file="regucalcin.pdf")
mygene1 <- getGene(cuff_data,'regucalcin')
expressionBarplot(mygene1)
dev.off()

expressionBarplot(isoforms(mygene1))

...

 
 Code Block
title
To plot gene level and isoform level expression for gene Rala
pdf(file="rala.pdf")
mygene2 <- getGene(cuff_data, 'Rala')
expressionBarplot(mygene2)
expressionBarplot(isoforms(mygene2))

dev.off()
...