Mercurial > repos > sblanck > mpagenomics
diff segmentation.R @ 4:3fcbb8030fcc draft
"planemo upload for repository https://github.com/sblanck/MPAgenomics4Galaxy/tree/master/mpagenomics_wrappers commit 40eda5ea3551e8b3bae32d0a8f405fe90ef22646-dirty"
author | sblanck |
---|---|
date | Mon, 12 Apr 2021 14:47:09 +0000 |
parents | 4d539083cf7f |
children |
line wrap: on
line diff
--- a/segmentation.R Tue Jun 16 04:34:09 2020 -0400 +++ b/segmentation.R Mon Apr 12 14:47:09 2021 +0000 @@ -3,7 +3,7 @@ options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } ) # we need that to not crash galaxy with an UTF8 error on German LC settings. -loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8") +# loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8") library("optparse") @@ -58,9 +58,13 @@ #signalType=args[8] library(MPAgenomics) -workdir=file.path(tmp_dir,"mpagenomics",userId) +workdir=file.path(tmp_dir) +if (!dir.exists(workdir)) + dir.create(workdir, showWarnings = TRUE, recursive = TRUE) setwd(workdir) +workdir + if (outputlog){ sinklog <- file(log, open = "wt") sink(sinklog ,type = "output") @@ -92,6 +96,21 @@ callobj= callingObject(copynumber=currentSeg$signal, segmented=currentSeg$segmented,chromosome=rep(chr,length(currentSeg$signal)), position=currentSeg$startPos,sampleNames=sample) currentCall=callingProcess(callobj,nclass=nbcall,cellularity=cellularity,verbose=TRUE) currentResult=currentCall$segment + if(outputgraph) + { + + currentPos=unlist(currentPositions) + figName <- sprintf("%s,%s", sample, chr); + pathname <- file.path(sprintf("%s.png", figName)); + png(filename = pathname, width = 1280, height = 480) + plot(NA,xlim=c(min(currentPos),max(currentPos)), ylim=c(0,6),xlab="Position", main=figName,ylab="CN", pch=".") + points(currentPos, unlist(currentSignal), pch="."); + for(i in 1:nrow(currentResult)) + lines(c(currentResult$chromStart[i],currentResult$chromEnd[i]),rep(currentResult$means[i],2),col="red",lwd=3) + dev.off() + + } + currentResult["sampleNames"]=c(rep(sample,length(currentCall$segment$chrom))) result=rbind(result,currentResult) } @@ -119,17 +138,38 @@ currentResult["chrom"]=c(rep(chr,length(currentSeg$segment$means))) currentResult["sampleNames"]=c(rep(sample,length(currentSeg$segment$means))) result=rbind(result,currentResult) - + if(outputgraph) + { + + currentPos=unlist(currentPositions) + figName <- sprintf("%s,%s", sample, chr); + pathname <- file.path(sprintf("%s.png", figName)); + png(filename = pathname, width = 1280, height = 480) + plot(NA,xlim=c(min(currentPos),max(currentPos)), ylim=c(0,1),xlab="Position", main=figName,ylab="CN", pch=".") + points(currentPos, unlist(currentSignal), pch="."); + print(currentResult) + for(i in 1:nrow(currentResult)) + lines(c(currentResult$start[i],currentResult$end[i]),rep(currentResult$means[i],2),col="red",lwd=3) + dev.off() + + } + + + } cat(paste0("OK\n")) } } finalResult=data.frame(sampleNames=result["sampleNames"],chrom=result["chrom"],chromStart=result["start"],chromEnd=result["end"],probes=result["points"],means=result["means"],stringsAsFactors=FALSE) + colnames(finalResult)=c("sampleNames","chrom","chromStart","chromEnd","probes","means") write.table(finalResult,output,row.names = FALSE, quote=FALSE, sep = "\t") } if (outputgraph){ - file.rename(file.path(tmp_dir,"mpagenomics",userId,"Rplots.pdf"), graph) + library(zip) + files2zip <- dir(pattern=".png") + zipr(graph, files = files2zip) + } if (outputlog){