view test-data/expressionplot.txt @ 5:78fcfc04fcfe draft

planemo upload for repository https://github.com/galaxyproject/tools-devteam/tree/master/tools/cummerbund commit f3f0bef4a450aafab3c6b05a27647471f93b22f3
author devteam
date Thu, 23 Feb 2017 20:24:03 -0500
parents 587c425b4e76
children
line wrap: on
line source

> get_features <- function(myGenes, f="gene") {
+     if (f == "isoforms")
+         return(isoforms(myGenes))
+     else if (f == "tss")
+         return(TSS(myGenes))
+     else if (f == "cds")
+         return(CDS(myGenes))
+     else
+         return(myGenes)
+ }
> 
> ## Main Function ##
> 
> library(argparse)
Loading required package: proto
> 
> parser <- ArgumentParser(description='Create a plot with cummeRbund')
> 
> parser$add_argument('--type', dest='plotType', default='Density', required=TRUE)
> parser$add_argument('--height', dest='height', type='integer', default=960, required=TRUE)
> parser$add_argument('--width', dest='width', type='integer', default=1280, required=TRUE)
> parser$add_argument('--outfile', dest='filename', default="plot-unknown-0.png", required=TRUE)
> parser$add_argument('--input', dest='input_database', default="cuffData.db", required=TRUE)
> parser$add_argument('--smooth', dest='smooth', action="store_true", default=FALSE)
> parser$add_argument('--gene_selector', dest='gene_selector', action="store_true", default=FALSE)
> parser$add_argument('--replicates', dest='replicates', action="store_true", default=FALSE)
> parser$add_argument('--labcol', dest='labcol', action="store_true", default=FALSE)
> parser$add_argument('--labrow', dest='labrow', action="store_true", default=FALSE)
> parser$add_argument('--border', dest='border', action="store_true", default=FALSE)
> parser$add_argument('--summary', dest='summary', action="store_true", default=FALSE)
> parser$add_argument('--count', dest='count', action="store_true", default=FALSE)
> parser$add_argument('--error_bars', dest='error_bars', action="store_true", default=FALSE)
> parser$add_argument('--log10', dest='log10', action="store_true", default=FALSE)
> parser$add_argument('--features', dest='features', action="store", default="genes")
> parser$add_argument('--clustering', dest='clustering', action="store", default="both")
> parser$add_argument('--iter_max', dest='iter_max', action="store")
> parser$add_argument('--genes', dest='genes', action="append")
> parser$add_argument('--k', dest='k', action="store")
> parser$add_argument('--x', dest='x', action="store")
> parser$add_argument('--y', dest='y', action="store")
> 
> args <- parser$parse_args()
> 
> ## Load cummeRbund library
> library("cummeRbund")
Loading required package: BiocGenerics
Loading required package: methods
Loading required package: parallel

Attaching package: ‘BiocGenerics’

The following objects are masked from ‘package:parallel’:

    clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
    clusterExport, clusterMap, parApply, parCapply, parLapply,
    parLapplyLB, parRapply, parSapply, parSapplyLB

The following objects are masked from ‘package:stats’:

    IQR, mad, xtabs

The following objects are masked from ‘package:base’:

    anyDuplicated, append, as.data.frame, cbind, colnames, do.call,
    duplicated, eval, evalq, Filter, Find, get, grep, grepl, intersect,
    is.unsorted, lapply, lengths, Map, mapply, match, mget, order,
    paste, pmax, pmax.int, pmin, pmin.int, Position, rank, rbind,
    Reduce, rownames, sapply, setdiff, sort, table, tapply, union,
    unique, unsplit, which, which.max, which.min

Loading required package: RSQLite
Loading required package: DBI
Loading required package: ggplot2
Loading required package: reshape2
Loading required package: fastcluster

Attaching package: ‘fastcluster’

The following object is masked from ‘package:stats’:

    hclust

Loading required package: rtracklayer
Loading required package: GenomicRanges
Loading required package: stats4
Loading required package: S4Vectors

Attaching package: ‘S4Vectors’

The following objects are masked from ‘package:base’:

    colMeans, colSums, expand.grid, rowMeans, rowSums

Loading required package: IRanges
Loading required package: GenomeInfoDb
Loading required package: Gviz
Loading required package: grid
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘ensembldb’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’
Warning: namespace ‘AnnotationHub’ is not available and has been replaced
by .GlobalEnv when processing object ‘plot.index’

Attaching package: ‘cummeRbund’

The following object is masked from ‘package:GenomicRanges’:

    promoters

The following object is masked from ‘package:IRanges’:

    promoters

The following object is masked from ‘package:BiocGenerics’:

    conditions

> 
> ## Initialize cuff object
> cuff <- readCufflinks(dir = "", dbFile = args$input_database, rebuild = FALSE)
> 
> ## Print out info
> print(cuff)
CuffSet instance with:
	 2 samples
	 87 genes
	 90 isoforms
	 88 TSS
	 0 CDS
	 87 promoters
	 88 splicing
	 0 relCDS
> sink("cuffdb_info.txt")
> print(cuff)
> print("SAMPLES:")
> samples(cuff)
> print("REPLICATES:")
> replicates(cuff)
> print("FEATURES:")
> print(annotation(genes(cuff)))
> cat(annotation(genes(cuff))[[1]],sep=",")
> sink()
> 
> png(filename = args$filename, width = args$width, height = args$height, type=c('cairo-png'))
> tryCatch({
+     if (args$plotType == 'density') {
+         csDensity(genes(cuff), replicates=args$replicates, logMode=args$log10)
+     }
+     else if (args$plotType == 'boxplot') {
+         csBoxplot(genes(cuff), replicates=args$replicates, logMode=args$log10)
+     }
+     else if (args$plotType == 'mds') {
+         MDSplot(genes(cuff), replicates=args$replicates)
+     }
+     else if (args$plotType == 'pca') {
+         PCAplot(genes(cuff), "PC1", "PC2", replicates=args$replicates)
+     }
+     else if (args$plotType == 'dendrogram') {
+         csDendro(genes(cuff), replicates=args$replicates)
+     }
+     else if (args$plotType == 'scatter') {
+         if (args$gene_selector) {
+             myGenes <- getGenes(cuff, args$genes)
+             csScatter(get_features(myGenes, args$features), args$x, args$y, smooth=args$smooth, logMode=args$log10)
+         }
+         else {
+             csScatter(genes(cuff), args$x, args$y, smooth=args$smooth, logMode=args$log10)
+         }
+     }
+     else if (args$plotType == 'volcano') {
+         if (args$gene_selector) {
+             myGenes <- get_features(getGenes(cuff, args$genes), args$features)
+         }
+         else {
+             myGenes <- genes(cuff)
+         }
+         csVolcano(myGenes, args$x, args$y)
+     }
+     else if (args$plotType == 'heatmap') {
+         if (args$gene_selector) {
+             myGenes <- getGenes(cuff, args$genes)
+         }
+         else {
+             myGenes <- getGenes(cuff,annotation(genes(cuff))[[1]])
+         }
+         csHeatmap(get_features(myGenes, args$features), clustering=args$clustering, labCol=args$labcol, labRow=args$labrow, border=args$border, logMode=args$log10)
+     }
+     else if (args$plotType == 'cluster') {
+         myGenes <- getGenes(cuff, args$genes)
+         csCluster(get_features(myGenes, args$features), k=args$k)
+     }
+     else if (args$plotType == 'dispersion') {
+         dispersionPlot(genes(cuff))
+     }
+     else if (args$plotType == 'fpkmSCV') {
+         fpkmSCVPlot(genes(cuff))
+     }
+     else if (args$plotType == 'scatterMatrix') {
+         csScatterMatrix(genes(cuff))
+     }
+     else if (args$plotType == 'expressionplot') {
+         myGenes <- getGenes(cuff, args$genes)
+         expressionPlot(get_features(myGenes, args$features), drawSummary=args$summary, showErrorbars=args$error_bars, replicates=args$replicates)
+     }
+     else if (args$plotType == 'expressionbarplot') {
+         myGeneId <- args$genes
+         myGenes <- getGenes(cuff, myGeneId)
+         expressionBarplot(get_features(myGenes, args$features), showErrorbars=args$error_bars, replicates=args$replicates)
+     }
+     else if (args$plotType == 'mds') {
+         MDSplot(genes(cuff),replicates=args$replicates)
+     }
+     else if (args$plotType == 'pca') {
+         PCAplot(genes(cuff),"PC1","PC2", replicates=args$replicates)
+     }
+     else if (args$plotType == 'maplot') {
+         MAplot(genes(cuff), args$x, args$y, useCount=args$count)
+     }
+     else if (args$plotType == 'genetrack') {
+         myGene <- getGene(cuff, args$genes)
+         plotTracks(makeGeneRegionTrack(myGene))
+     }
+ },error = function(e) {
+     write(paste("Failed:", e, sep=" "), stderr())
+     q("no", 1, TRUE)
+ })
Getting gene information:
	FPKM
	Differential Expression Data
	Annotation Data
	Replicate FPKMs
	Counts
Getting isoforms information:
	FPKM
	Differential Expression Data
	Annotation Data
	Replicate FPKMs
	Counts
Getting CDS information:
	FPKM
	Differential Expression Data
	Annotation Data
	Replicate FPKMs
	Counts
Getting TSS information:
	FPKM
	Differential Expression Data
	Annotation Data
	Replicate FPKMs
	Counts
Getting promoter information:
	distData
Getting splicing information:
	distData
Getting relCDS information:
	distData
> devname = dev.off()
>