User Tools

Site Tools


mkatari-bioinformatics-august-2013-deseq

Differences

This shows you the differences between two versions of the page.

Link to this comparison view

Both sides previous revisionPrevious revision
Next revision
Previous revision
mkatari-bioinformatics-august-2013-deseq [2015/06/17 06:15] mkatarimkatari-bioinformatics-august-2013-deseq [2015/08/21 14:13] (current) mkatari
Line 21: Line 21:
 First we will load the count data file. First we will load the count data file.
 <code> <code>
-counts = read.table(pathToCountsData, header=T, row.names=1)+counts = read.table("NextGenRaw.txt", header=T, row.names=1)
 </code> </code>
  
 Then we will load the experimental design. An example is provided [[https://docs.google.com/file/d/0B172nc4dAaaOaE5fTVVhUHJKazg/edit?usp=sharing|here]]: Then we will load the experimental design. An example is provided [[https://docs.google.com/file/d/0B172nc4dAaaOaE5fTVVhUHJKazg/edit?usp=sharing|here]]:
 <code> <code>
-expdesign = read.table(pathToExpDesign)+expdesign = read.table("expdesign.txt")
 </code> </code>
  
Line 53: Line 53:
 An important part of DESeq is to estimate dispersion. This is simply a form of variance for the genes. An important part of DESeq is to estimate dispersion. This is simply a form of variance for the genes.
 <code> <code>
 +# if you have replicates do the following:
 cds = estimateDispersions( cds ) cds = estimateDispersions( cds )
 +### HOWEVER, If you have NO replicates, then try this
 +cds = estimateDispersions( cds, method="blind" , sharingMode = "fit-only" )
 </code> </code>
  
 To visualize the disperson graph To visualize the disperson graph
 <code> <code>
-dispersionFile = paste(pathToOutputDir, "Dispersion.pdf", sep="")+dispersionFile = "Dispersion.pdf"
 pdf(dispersionFile) pdf(dispersionFile)
 plotDispEsts( cds ) plotDispEsts( cds )
Line 77: Line 80:
  
 <code> <code>
-maFile = paste(pathToOutputDir, "MAplot.pdf", sep="")+maFile = "MAplot.pdf"
 pdf(maFile) pdf(maFile)
 plotMA(res) plotMA(res)
Line 85: Line 88:
 #To get the genes that have FDR of 10% and save it in the output directory. #To get the genes that have FDR of 10% and save it in the output directory.
 <code> <code>
-resSig = res[ which(res$padj < 0.1), ]+resSigind = res[ which(res$padj < 0.1 & res$log2FoldChange > 1), ] 
 +resSigrep = res[ which(res$padj < 0.1 & res$log2FoldChange < -1), ]
  
-outfile paste(pathToOutputDir,"Deseq.results.txt", sep=""+indoutfile = "Deseq.indresults.txt"  
 +repoutfile = "Deseq.represults.txt
  
-write.table(resSig,  +write.table(resSigind,  
-            outfile+            indoutfile,  
 +            sep="\t",  
 +            col.names=T,  
 +            row.names=F, 
 +            quote=F) 
 + 
 +write.table(resSigrep,  
 +            repoutfile
             sep="\t",              sep="\t", 
             col.names=T,              col.names=T, 
mkatari-bioinformatics-august-2013-deseq.1434521738.txt.gz · Last modified: 2015/06/17 06:15 by mkatari