%\VignetteIndexEntry{DEGreport} %\VignetteKeywords{DifferentialExpression, Visualization, RNASeq, ReportWriting} %\VignetteEngine{knitr::knitr} \documentclass[12pt]{article} \newcommand{\deseqtwo}{\textit{DESeq2}} \newcommand{\lowtilde}{\raise.17ex\hbox{$\scriptstyle\mathtt{\sim}$}} <>= library("knitr") opts_chunk$set(tidy=FALSE,dev="png",fig.show="hide", fig.width=4,fig.height=4.5, message=FALSE) @ <>= BiocStyle::latex() @ \title{DEGreport } \author{Lorena Pantano} \date{Modified: 23 May, 2014. Compiled: \today} \begin{document} \maketitle <>= library(DEGreport) data(humanSexDEedgeR) library(edgeR) @ We are going to do a differential expression analysis with edgeR. We have an object that is comming from the edgeR package. It countains a gene count matrix for 85 TSI HapMap individuals, and the gender information. With that, we are going to apply the `glmFit` function to get genes differentially expressed between males and females. <>= des<-humanSexDEedgeR$design fit <- glmFit(humanSexDEedgeR,des) lrt <- glmLRT(fit) tab<-cbind(lrt$table,p.adjust(lrt$table$PValue,method="BH")) detags <- rownames(tab[tab[,5]<=0.1,]) plotSmear(humanSexDEedgeR, de.tags=detags) @ We need to extract the experiment design data.frame where the condition is Male or Female. <>= counts<-cpm(humanSexDEedgeR,log=FALSE) g1<-colnames(counts)[1:41] g2<-colnames(counts)[42:85] design<-data.frame(condition=sub("1","Male",sub("0","Female",des[,2]))) @ We are getting the chromosome information for each gene. This way we can colour genes according autosomic,X or Y chromosomes. <>= data(geneInfo) @ Create the report. The main parameters are the column names in group1, and group2. Then, the count matrix, gene names that are DE, p-values, fold changes and path to create the report. As optional, you can give colours for each gene, and the number of permutation. <>= detag10<-detags[1:10] pval<-tab[,4] fc<-tab[detag10,1] @ Run the following lines to create the file <>= pathreport<-"~/report" #change this to a proper path createReport(g1,g2,counts,detag10,pval,fc,pathreport,colors,pop=400) @ Run the fowlling lines if you want to visualize your expression values by condition: <>= degObj(counts,design,"/tmp/degObj.rda") library(shiny) runGist(9930881) @ You can use individual functions, like degRank or degMean. This will create specific figures and tables that are included in the report. <>= degMean(pval,counts) degVar(pval,counts) degMV(g1,g2,pval,counts) degMB(detags,g1,g2,counts) degVB(detags,g1,g2,counts) @ <>== library(coda) library(rjags) rank<-degRank(g1,g2,counts[detag10,],fc,400,500) degPR(rank) @ \end{document}