Mercurial > repos > stef > qdnaseq_test

comparison QDNAseq-plot.R @ 27:f89205f51e27 draft default tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author stef
date Mon, 06 Jul 2015 05:41:08 -0400
parents 4943308e95fc
children
comparison
equal deleted inserted replaced
26:069289631381 27:f89205f51e27
35 rVersion <- R.version.string 35 rVersion <- R.version.string
36 catMsg( c("QDNAseq version: ", qdnaseqVersion) ) 36 catMsg( c("QDNAseq version: ", qdnaseqVersion) )
37 catMsg( c( rVersion ) ) 37 catMsg( c( rVersion ) )
38 38
39 qdnaseqObject <- readRDS( rdsFilePath ) 39 qdnaseqObject <- readRDS( rdsFilePath )
40 chromosomesToPlot <- unlist( strsplit( chromosomesToPlotString, ",") )
41
42 #cat( "CHROM: ", chromosomesToPlotString, "\n" )
43 #cat( "REGION: ", regionToPlotString, "\n" )
44 #cat( "What to plot: ", whatToPlot, "\n" )
45
40 ## COPYNUMBER PLOT 46 ## COPYNUMBER PLOT
41 sample <- SAMPLE_INDEX 47 sample <- SAMPLE_INDEX
42 png( outputPngPath, width=PLOT_WIDTH, height=PLOT_HEIGHT, res=PLOT_RES ) 48 png( outputPngPath, width=PLOT_WIDTH, height=PLOT_HEIGHT, res=PLOT_RES )
43 par( PAR_SET ) 49 par( PAR_SET )
44 plot( qdnaseqObject[ ,sample ] ) 50 ylab_text <- "log2 read counts"
51 if ( whatToPlot == 'everything' ){
52 catMsg( c( "Plotting all data in object" ) )
53 plot( qdnaseqObject[ ,sample ], ylab=ylab_text )
54 abline( h=c(-4,-3,-2,-1,1,2,3,4), lty=1, lwd=0.5, col="grey" )
55 } else if( whatToPlot == 'chromosomes' ){
56 catMsg( c( "Plotting subset of chromosomes" ) )
57 fdata <- qdnaseqObject@featureData@data
58 idx_region <- which( fdata$chromosome %in% chromosomesToPlot )
59 plot( qdnaseqObject[ idx_region, sample ], ylab=ylab_text )
60 abline( h=c(-4,-3,-2,-1,1,2,3,4), lty=1, lwd=0.5, col="grey" )
61 } else if( whatToPlot == 'region' ){
62 regionC <- chrName
63 regionS <- chrStart
64 regionE <- chrEnd
65 if ( regionS > regionE ) stop("Chosen start is > end")
66 catMsg( c( "Plotting genomic region (chr=", regionC, " start=", regionS, " end=", regionE, ")" ) )
67
68 fdata <- qdnaseqObject@featureData@data
69 idx_region <- which( fdata$chromosome == regionC & fdata$start > regionS & fdata$end < regionE )
70 cat( idx_region, "\n")
71
72 plot( qdnaseqObject[ idx_region, sample ], doCalls=FALSE, ylab=ylab_text )
73 }
45 #mtext( "plotted in galaxy", 3 ) 74 #mtext( "plotted in galaxy", 3 )
46 abline( h=c(-2,-1,1,2,3,4), lty=1, lwd=0.5, col="grey" ) 75
47 dev.off() 76 dev.off()
48 77
49 78
50 ## all ok 79 ## all ok
51 q(status=0) 80 q(status=0)