Mercurial > repos > fabio > iwtomics

comparison loadandplot.R @ 83:f81d72e482cf draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author fabio
date Fri, 06 Oct 2017 11:06:44 -0400
parents
children
comparison
equal deleted inserted replaced
82:8e397b62df4c 83:f81d72e482cf
1 if (require("IWTomics",character.only = TRUE,quietly = FALSE)) {
2 args=commandArgs(TRUE)
3
4 # get args names and values
5 args_values=strsplit(args,'=')
6 args_names=unlist(lapply(args_values,function(arg) arg[1]))
7 names(args_values)=args_names
8 args_values=lapply(args_values,function(arg) arg[2])
9 # read filenames
10 outrdata=args_values$outrdata
11 outregions=args_values$outregions
12 outfeatures=args_values$outfeatures
13 outpdf=args_values$outpdf
14 regionspaths=unlist(strsplit(args_values$regionspaths,'\\|'))
15 if("regionsheaderfile" %in% args_names){
16 # the file regionsheaderfile must contain as first column the (unique) regionsfilenames,
17 # as second column the corresponding ids and as third column the names
18 tryCatch({
19 regionsheader=read.delim(args_values$regionsheaderfile,header=FALSE,stringsAsFactors=FALSE,row.names=1,sep="\t")
20 regionsfilenames=unlist(strsplit(args_values$regionsfilenames,'\\|'))
21 if(length(setdiff(regionsfilenames,row.names(regionsheader)))) {
22 write("IWTomics message: Not all region files are present in the first column of header file for regions.", stderr())
23 quit(save="no", status=11)
24 }
25 id_regions=regionsheader[regionsfilenames,1]
26 name_regions=regionsheader[regionsfilenames,2]
27 }, error = function(err) {
28 write("IWTomics message: An error has occurred reading the header file for regions. Please try again.", stderr())
29 quit(save="no", status=10) #error on header file
30 })
31 }else{
32 eval(parse(text=args[[which(args_names=='regionsgalaxyids')]]))
33 id_regions=paste0('data_',regionsgalaxyids)
34 name_regions=paste0('data_',regionsgalaxyids)
35 }
36 featurespaths=unlist(strsplit(args_values$featurespaths,'\\|'))
37 if("featuresheaderfile" %in% args_names){
38 # the file featuresheaderfile must contain as first column the (unique) featuresfilenames,
39 # as second column the corresponding ids and as third column the names
40 tryCatch({
41 featuresheader=read.delim(args_values$featuresheaderfile,header=FALSE,stringsAsFactors=FALSE,row.names=1,sep="\t")
42 featuresfilenames=unlist(strsplit(args_values$featuresfilenames,'\\|'))
43 if(length(setdiff(featuresfilenames,row.names(featuresheader)))) {
44 write("IWTomics message: Not all feature files are present in the first column of header file for features.", stderr())
45 quit(save="no", status=21)
46 }
47 id_features=featuresheader[featuresfilenames,1]
48 name_features=featuresheader[featuresfilenames,2]
49 }, error = function(err) {
50 write("IWTomics message: An error has occurred reading the header file for features. Please try again.", stderr())
51 quit(save="no", status=20) #error on header file
52 })
53 }else{
54 eval(parse(text=args[[which(args_names=='featuresgalaxyids')]]))
55 id_features=paste0('data_',featuresgalaxyids)
56 name_features=paste0('data_',featuresgalaxyids)
57 }
58 # read parameters (from smoothing on)
59 i_smoothing=which(args_names=='smoothing')
60 for(i in i_smoothing:length(args)){
61 eval(parse(text=args[[i]]))
62 }
63
64 # load data
65 tryCatch({
66 regionsFeatures=IWTomicsData(regionspaths,featurespaths,alignment,
67 id_regions,name_regions,id_features,name_features,start.are.0based=start.are.0based)
68 }, error = function(err) {
69 if(grepl('invalid format',err$message)){
70 write("IWTomics message: Not enough columns in input file.", stderr())
71 quit(save="no", status=31) # error, not enough columns in input file
72 }else if(grepl('duplicated regions',err$message)){
73 write("IWTomics message: Duplicated regions in region file.", stderr())
74 quit(save="no", status=32) # error, duplicated regions in region file
75 }else if(grepl('duplicated windows',err$message)){
76 write("IWTomics message: Duplicated windows in feature file.", stderr())
77 quit(save="no", status=33) # error, duplicated windows in feature file
78 }else if(grepl('overlapping windows',err$message)){
79 write("IWTomics message: Overlapping windows in feature file.", stderr())
80 quit(save="no", status=34) # error, overlapping windows in feature file
81 }else if(grepl('not all regions in datasets',err$message)){
82 write("IWTomics message: Windows in feature files do not cover all regions in region files.", stderr())
83 quit(save="no", status=35) # error, windows in feature files do not cover all regions in region files
84 }else if(grepl('ifferent size windows',err$message)){
85 write("IWTomics message: All windows in a feature file must have the same size.", stderr())
86 quit(save="no", status=36) # error, all windows in a feature files must have the same size
87 }
88 #error loading data
89 write("IWTomics message: An error has occurred reading the data. Please try again.", stderr())
90 quit(save="no", status=30)
91 })
92
93 # smooth data
94 if(smoothing!='no'){
95 tryCatch({
96 if(smoothing=='locpoly'){
97 dist_knots=10
98 }else if(smoothing=='kernel'){
99 degree=3
100 dist_knots=10
101 }else if(smoothing=='splines'){
102 bandwidth=5
103 }
104 if(alignment=='scale'){
105 if(scale==0){
106 regionsFeatures=smooth(regionsFeatures,type=smoothing,fill_gaps=fill_gaps,
107 bandwidth=bandwidth,degree=degree,dist_knots=dist_knots)
108 }else{
109 regionsFeatures=smooth(regionsFeatures,type=smoothing,fill_gaps=fill_gaps,
110 bandwidth=bandwidth,degree=degree,dist_knots=dist_knots,scale_grid=scale)
111 }
112 }else{
113 regionsFeatures=smooth(regionsFeatures,type=smoothing,fill_gaps=fill_gaps,
114 bandwidth=bandwidth,degree=degree,dist_knots=dist_knots)
115 }
116 }, error = function(err) {
117 write("IWTomics message: An error has occurred smoothing the data. Please try again.", stderr())
118 quit(save="no", status=40) #error on smoothing
119 })
120 }
121
122 # plot data
123 pdf(outpdf,width=10,height=8)
124 if(plottype=='boxplot'){
125 # fix repeated probs
126 probs=sort(unique(probs))
127 }else{
128 probs=c(0.25,0.5,0.75)
129 }
130 plot(regionsFeatures,type=plottype,probs=probs,average=average,size=size,ask=FALSE)
131 dev.off()
132
133 # create output
134 #write.table(cbind(unlist(strsplit(args_values$regionsfilenames,'\\|')),idRegions(regionsFeatures),nameRegions(regionsFeatures)),
135 #file=outregions,quote=FALSE,sep='\t',row.names=FALSE,col.names=FALSE)
136 write.table(as.data.frame(t(idRegions(regionsFeatures))),file=outregions,quote=FALSE,sep='\t',row.names=FALSE,col.names=FALSE)
137 #write.table(cbind(unlist(strsplit(args_values$featuresfilenames,'\\|')),idFeatures(regionsFeatures),nameFeatures(regionsFeatures)),
138 #file=outfeatures,quote=FALSE,sep='\t',row.names=FALSE,col.names=FALSE)
139 write.table(as.data.frame(t(idFeatures(regionsFeatures))),file=outfeatures,quote=FALSE,sep='\t',row.names=FALSE,col.names=FALSE)
140 save(regionsFeatures,file=outrdata)
141 }else{
142 write("IWTomics message: Missing IWTomics package. Please be sure to have it installed before using this tool.", stderr())
143 quit(save="no", status=255)
144 }