Mercurial > repos > sblanck > mpagenomics

comparison extractCN.R @ 4:3fcbb8030fcc draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"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
comparison
equal deleted inserted replaced
3:94fc6ed13946 4:3fcbb8030fcc
1 #!/usr/bin/env Rscript 1 #!/usr/bin/env Rscript
2 # setup R error handling to go to stderr 2 # setup R error handling to go to stderr
3 options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } ) 3 options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } )
4 4
5 # we need that to not crash galaxy with an UTF8 error on German LC settings. 5 # we need that to not crash galaxy with an UTF8 error on German LC settings.
6 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8") 6 # loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
7 7
8 library("optparse") 8 library("optparse")
9 library("zip")
9 10
10 ##### Read options 11 ##### Read options
11 option_list=list( 12 option_list=list(
12 make_option("--chrom",type="character",default=NULL, dest="chrom"), 13 make_option("--chrom",type="character",default=NULL, dest="chrom"),
13 make_option("--input",type="character",default=NULL, dest="input"), 14 make_option("--input",type="character",default=NULL, dest="input"),
15 make_option("--zip",type="character",default=NULL, dest="zip"),
14 make_option("--output",type="character",default=NULL, dest="output"), 16 make_option("--output",type="character",default=NULL, dest="output"),
15 make_option("--new_file_path",type="character",default=NULL, dest="new_file_path"), 17 make_option("--new_file_path",type="character",default=NULL, dest="new_file_path"),
16 make_option("--settings_type",type="character",default=NULL, dest="settings_type"), 18 make_option("--settings_type",type="character",default=NULL, dest="settings_type"),
17 make_option("--settings_tumor",type="character",default=NULL, dest="settings_tumor"), 19 make_option("--settings_tumor",type="character",default=NULL, dest="settings_tumor"),
18 make_option("--symmetrize",type="character",default=NULL, dest="symmetrize"), 20 make_option("--symmetrize",type="character",default=NULL, dest="symmetrize"),
33 35
34 #loading libraries 36 #loading libraries
35 37
36 chrom=opt$chrom 38 chrom=opt$chrom
37 input=opt$input 39 input=opt$input
40 zip=opt$zip
38 tmp_dir=opt$new_file_path 41 tmp_dir=opt$new_file_path
39 output=opt$output 42 output=opt$output
40 settingsType=opt$settings_type 43 settingsType=opt$settings_type
41 tumorcsv=opt$settings_tumor 44 tumorcsv=opt$settings_tumor
42 symmetrize=opt$symmetrize 45 symmetrize=opt$symmetrize
45 outputlog=opt$outputlog 48 outputlog=opt$outputlog
46 log=opt$log 49 log=opt$log
47 user=opt$userid 50 user=opt$userid
48 51
49 library(MPAgenomics) 52 library(MPAgenomics)
50 workdir=file.path(tmp_dir, "mpagenomics",user) 53 library(aroma.affymetrix)
54 library(R.utils)
55 #workdir=file.path(tmp_dir, "mpagenomics",user)
56 tmp_dir
57 tmp_dir=file.path(tmp_dir)
58 if (!dir.exists(tmp_dir))
59 dir.create(tmp_dir, showWarnings = TRUE, recursive = TRUE)
60 setwd(tmp_dir)
61 # tmpzip=file.copy(from = zip,to=paste0(workdir,"/tmp.zip"))
62 # tmpzip
63 unzip(zipfile = zip,exdir = ".")
64 # if (file.exists(tmpzip)) {
65 # #Delete file if it exists
66 # file.remove(fn)
67 # }
68
69 workdir=file.path(tmp_dir,user)
51 setwd(workdir) 70 setwd(workdir)
52
53 inputDataset=read.table(file=input,stringsAsFactors=FALSE) 71 inputDataset=read.table(file=input,stringsAsFactors=FALSE)
54 dataset=inputDataset[1,2] 72 dataset=inputDataset[1,2]
55 73
56 if (outputlog){ 74 if (outputlog){
57 sinklog <- file(log, open = "wt") 75 sinklog <- file(log, open = "wt")
73 if (tumorcsv== "None") 91 if (tumorcsv== "None")
74 { 92 {
75 CN=getCopyNumberSignal(dataset,chromosome=chrom_vec, onlySNP=snp) 93 CN=getCopyNumberSignal(dataset,chromosome=chrom_vec, onlySNP=snp)
76 94
77 } else { 95 } else {
78 CN=getCopyNumberSignal(dataset,chromosome=chrom_vec, normalTumorArray=tumorcsv, onlySNP=snp) 96 CN=getCopyNumberSignal(dataset,chromosome=chrom_vec, normalTumorArray=tumorcsv, onlySNP=snp)
79 } 97 }
80 } else { 98 } else {
81 input_tmp <- strsplit(settingsType,",") 99 input_tmp <- strsplit(settingsType,",")
82 input_tmp_vecstring <-unlist(input_tmp) 100 input_tmp_vecstring <-unlist(input_tmp)
83 input_vecstring = sub("^([^.]*).*", "\\1", input_tmp_vecstring) 101 input_vecstring = sub("^([^.]*).*", "\\1", input_tmp_vecstring)
107 input_tmp_vecstring <-unlist(input_tmp) 125 input_tmp_vecstring <-unlist(input_tmp)
108 input_vecstring = sub("^([^.]*).*", "\\1", input_tmp_vecstring) 126 input_vecstring = sub("^([^.]*).*", "\\1", input_tmp_vecstring)
109 } 127 }
110 128
111 symFracB_global=data.frame(check.names = FALSE) 129 symFracB_global=data.frame(check.names = FALSE)
112 130 tumorFile=read.csv(tumorcsv,header=TRUE)
131 tumor=tumorFile$tumor
132 input_vecstring=input_vecstring[which(input_vecstring %in% tumor)]
113 for (currentFile in input_vecstring) { 133 for (currentFile in input_vecstring) {
114 cat(paste0("extracting signal from ",currentFile,".\n")) 134 cat(paste0("extracting signal from ",currentFile,".\n"))
115 currentSymFracB=data.frame() 135 currentSymFracB=data.frame()
116 symFracB=getSymFracBSignal(dataset,chromosome=chrom_vec,file=currentFile,normalTumorArray=tumorcsv) 136 symFracB=getSymFracBSignal(dataset,chromosome=chrom_vec,file=currentFile,normalTumorArray=tumorcsv)
117 list_chr=names(symFracB) 137 list_chr=names(symFracB)
122 142
123 } 143 }
124 if (is.null(symFracB_global) || nrow(symFracB_global)==0) { 144 if (is.null(symFracB_global) || nrow(symFracB_global)==0) {
125 symFracB_global=currentSymFracB 145 symFracB_global=currentSymFracB
126 } else { 146 } else {
127 symFracB_global=cbind(symFracB_global,currentFile=currentSymFracB[[3]]) 147 #symFracB_global=cbind(symFracB_global,currentFile=currentSymFracB[[3]])
148
149 symFracB_global=merge(symFracB_global,currentSymFracB[,c(3,4)],by="featureNames")
150 symFracB_global=symFracB_global[c(2:ncol(symFracB_global),1)]
151 symFracB_global=symFracB_global[order(symFracB_global$chromosome, symFracB_global$position),]
128 } 152 }
129 } 153 }
130 names(symFracB_global)[names(symFracB_global)=="featureNames"] <- "probeName" 154 names(symFracB_global)[names(symFracB_global)=="featureNames"] <- "probeName"
131 155
132 write.table(format(symFracB_global), output, row.names = FALSE, quote = FALSE, sep = "\t") 156 write.table(format(symFracB_global), output, row.names = FALSE, quote = FALSE, sep = "\t")
161 write.table(format(fracB_global), output, row.names = FALSE, quote = FALSE, sep = "\t") 185 write.table(format(fracB_global), output, row.names = FALSE, quote = FALSE, sep = "\t")
162 } 186 }
163 187
164 } 188 }
165 189
190 if (dir.exists(workdir))
191 system(paste0("rm -r ", workdir))
192
166 if (outputlog){ 193 if (outputlog){
167 sink(type="output") 194 sink(type="output")
168 sink(type="message") 195 sink(type="message")
169 close(sinklog) 196 close(sinklog)
170 } 197 }