|
0
|
1 #!/usr/bin/env Rscript
|
|
|
2
|
|
5
|
3 suppressPackageStartupMessages(library("data.table"))
|
|
0
|
4 suppressPackageStartupMessages(library("optparse"))
|
|
|
5
|
|
|
6 option_list <- list(
|
|
61
|
7 make_option(c("--build"), action="store", dest="build", help="Genome build"),
|
|
|
8 make_option(c("--chrom_len_file"), action="store", dest="chrom_len_file", help="Chromosome length file"),
|
|
|
9 make_option(c("--email"), action="store", dest="email", help="User email address"),
|
|
|
10 make_option(c("--galaxy_url"), action="store", dest="galaxy_url", help="Galaxy instance base URL"),
|
|
77
|
11 make_option(c("--hub_long_label"), action="store", dest="hub_long_label", help="Hub long label"),
|
|
61
|
12 make_option(c("--hub_name"), action="store", dest="hub_name", default=NULL, help="Hub name without spaces"),
|
|
77
|
13 make_option(c("--hub_short_label"), action="store", dest="hub_short_label", help="Hub short label"),
|
|
66
|
14 make_option(c("--input_dir_para"), action="store", dest="input_dir_para", help="Directory containing .para outputs from IDEAS"),
|
|
61
|
15 make_option(c("--input_dir_state"), action="store", dest="input_dir_state", help="Directory containing .state outputs from IDEAS"),
|
|
|
16 make_option(c("--output_trackhub"), action="store", dest="output_trackhub", help="Output hub file"),
|
|
|
17 make_option(c("--output_trackhub_files_path"), action="store", dest="output_trackhub_files_path", help="Output hub extra files path"),
|
|
|
18 make_option(c("--output_trackhub_id"), action="store", dest="output_trackhub_id", help="Encoded output_trackhub dataset id"),
|
|
68
|
19 make_option(c("--script_dir"), action="store", dest="script_dir", help="R script source directory"),
|
|
77
|
20 make_option(c("--state_colors"), action="store", dest="state_colors", default=NULL, help="List of state colors"),
|
|
|
21 make_option(c("--state_indexes"), action="store", dest="state_indexes", default=NULL, help="List of state indexes"),
|
|
|
22 make_option(c("--state_names"), action="store", dest="state_names", default=NULL, help="List of state names")
|
|
0
|
23 )
|
|
|
24
|
|
|
25 parser <- OptionParser(usage="%prog [options] file", option_list=option_list)
|
|
|
26 args <- parse_args(parser, positional_arguments=TRUE)
|
|
|
27 opt <- args$options
|
|
|
28
|
|
78
|
29 create_primary_html = function(output_trackhub, tracks_dir, build) {
|
|
|
30 track_files <- list.files(path=tracks_dir);
|
|
54
|
31 s <- paste('<html><head></head><body>', sep="\n");
|
|
78
|
32 s <- paste(s, '<h3>Contents of directory ~/myHub/', build, ' required by UCSC TrackHub</h3>\n', sep="");
|
|
54
|
33 s <- paste(s, '<ul>\n', sep="")
|
|
78
|
34 for (i in 1:length(track_files)) {
|
|
|
35 s <- paste(s, '<li><a href="', 'myHub/', build, "/", track_files[i], '">', track_files[i], '</a></li>\n', sep="");
|
|
54
|
36 }
|
|
|
37 s <- paste(s, '</ul>\n</body>\n</html>', sep="");
|
|
|
38 cat(s, file=output_trackhub);
|
|
20
|
39 }
|
|
|
40
|
|
77
|
41 create_cells = function(input_dir_state, chrom_len_file, base_track_file_name, state_indexes, state_names, state_colors) {
|
|
66
|
42 # Create everything needed, including the bigbed file,
|
|
54
|
43 # to render the tracks within the UCSC track hub.
|
|
|
44 state_files <- list.files(path=input_dir_state, full.names=TRUE);
|
|
|
45 genome_size = read.table(chrom_len_file);
|
|
|
46 g = NULL;
|
|
|
47 for(i in state_files) {
|
|
|
48 tg = as.matrix(fread(i));
|
|
|
49 t = NULL;
|
|
|
50 for(j in 1:dim(genome_size)[1]) {
|
|
|
51 t = c(t, which(tg[,2]==as.character(genome_size[j, 1]) & as.numeric(tg[,4]) > as.numeric(genome_size[j, 2])));
|
|
|
52 }
|
|
|
53 if (length(t) > 0) {
|
|
|
54 tg = tg[-t,];
|
|
|
55 }
|
|
|
56 t = which(is.na(match(tg[,2], genome_size[,1]))==T);
|
|
|
57 if (length(t)>0) {
|
|
|
58 tg = tg[-t,];
|
|
|
59 }
|
|
|
60 g = rbind(g, tg);
|
|
|
61 }
|
|
|
62 uchr = sort(unique(as.character(g[,2])));
|
|
|
63 g1 = NULL;
|
|
|
64 for(i in uchr) {
|
|
|
65 t = which(g[,2]==i);
|
|
|
66 g1 = rbind(g1, g[t[order(as.integer(g[t, 3]))],]);
|
|
|
67 }
|
|
|
68 g = NULL;
|
|
|
69 chr = as.character(g1[,2]);
|
|
|
70 posst = as.numeric(g1[,3]);
|
|
|
71 posed = as.numeric(g1[,4]);
|
|
|
72 state = as.matrix(g1[,5:(dim(g1)[2]-1)]);
|
|
|
73 state_name = 0:max(state);
|
|
|
74 L = dim(g1)[1];
|
|
|
75 n = dim(state)[2];
|
|
|
76 cells = as.character(colnames(g1)[5:(dim(g1)[2]-1)]);
|
|
|
77 options(scipen=999);
|
|
|
78 tt = which(chr[2:L]!=chr[2:L-1]);
|
|
|
79 tt = c(tt, which(posst[2:L]!=posed[2:L-1]));
|
|
|
80 tt = sort(unique(tt));
|
|
|
81 for(i in 1:n) {
|
|
|
82 tstate = state[,i];
|
|
|
83 t = c(tt, which(tstate[2:L]!=tstate[2:L-1]));
|
|
|
84 t = sort(unique(t));
|
|
|
85 t0 = c(0, t) + 1;
|
|
|
86 t = c(t, L);
|
|
|
87 np = cbind(chr[t], posst[t0], posed[t], tstate[t]);
|
|
57
|
88 track_file_name_bed_unsorted <- get_track_file_name(base_track_file_name, i, "bed_unsorted");
|
|
58
|
89 track_file_name_bed <- get_track_file_name(base_track_file_name, i, "bed");
|
|
56
|
90 track_file_name_bigbed <- get_track_file_name(base_track_file_name, i, "bigbed");
|
|
54
|
91 x = cbind(np[, 1:3], state_name[as.integer(np[,4])+1], 1000, ".", np[,2:3]);
|
|
57
|
92 write.table(as.matrix(x), track_file_name_bed_unsorted, quote=F, row.names=F, col.names=F);
|
|
59
|
93 cmd = paste("LC_COLLATE=C sort -k1,1 -k2,2n < ", track_file_name_bed_unsorted, " > ", track_file_name_bed);
|
|
|
94 system(cmd);
|
|
|
95 cmd = paste("bedToBigBed ", track_file_name_bed, chrom_len_file, " ", track_file_name_bigbed);
|
|
|
96 system(cmd);
|
|
65
|
97 system(paste("rm ", track_file_name_bed_unsorted));
|
|
|
98 system(paste("rm ", track_file_name_bed));
|
|
54
|
99 }
|
|
|
100 return(cells);
|
|
0
|
101 }
|
|
|
102
|
|
77
|
103 create_track_db = function(galaxy_url, encoded_dataset_id, input_dir_para, input_dir_state, build, chrom_len_file,
|
|
|
104 tracks_dir, hub_name, hub_short_label, hub_long_label, state_indexes, state_names, state_colors) {
|
|
54
|
105 # Create a trackDb.txt file that includes each state.
|
|
80
|
106 s_indexes = c();
|
|
|
107 s_names = c();
|
|
|
108 s_colors = c();
|
|
|
109 if (is.null(state_indexes)) {
|
|
|
110 # Generate state colors automatically.
|
|
|
111 index = 0;
|
|
|
112 para_files = list.files(path=input_dir_para, full.names=TRUE);
|
|
|
113 for (para_file in para_files) {
|
|
|
114 data_frame = read.table(para_file, comment="!", header=T);
|
|
|
115 color_hex_code = create_heatmap(data_frame);
|
|
|
116 s_indexes = c(s_indexes, index);
|
|
|
117 s_names = c(s_names, toString(index));
|
|
|
118 s_colors = c(s_colors, color_hex_code);
|
|
|
119 }
|
|
|
120 } else {
|
|
61
|
121 # Split state_indexes into a list of integers.
|
|
80
|
122 index_str = as.character(state_indexes);
|
|
|
123 items = strsplit(index_str, ",")
|
|
61
|
124 for (item in items) {
|
|
80
|
125 s_indexes = c(s_indexes, as.integer(item));
|
|
61
|
126 }
|
|
77
|
127 # Split state_names into a list of strings.
|
|
80
|
128 name_str = as.character(state_colors);
|
|
|
129 items = strsplit(name_str, ",");
|
|
|
130 item_index = 0;
|
|
77
|
131 for (item in items) {
|
|
|
132 item_index = item_index + 1;
|
|
|
133 # Handle the special string "use state index".
|
|
|
134 if (identical(item, "use state index")) {
|
|
80
|
135 item = s_indexes[item_index];
|
|
77
|
136 }
|
|
80
|
137 s_names = c(s_names, item);
|
|
77
|
138 }
|
|
61
|
139 # Split state_colors into a list of strings.
|
|
80
|
140 color_str = as.character(state_colors);
|
|
|
141 items = strsplit(color_str, ",");
|
|
61
|
142 for (item in items) {
|
|
80
|
143 s_colors = c(s_colors, item);
|
|
61
|
144 }
|
|
|
145 }
|
|
77
|
146 base_track_file_name <- paste(tracks_dir, hub_name, sep="");
|
|
|
147 cells = create_cells(input_dir_state, chrom_len_file, base_track_file_name, s_indexes, s_names, s_colors);
|
|
54
|
148 track_db = NULL;
|
|
|
149 for (i in 1:length(cells)) {
|
|
71
|
150 # Get the color for the current state.
|
|
80
|
151 color <- paste(c(col2rgb(s_colors[i])), collapse=",");
|
|
71
|
152 # Get the bigDataUrl.
|
|
54
|
153 big_data_url <- get_big_data_url(galaxy_url, encoded_dataset_id, tracks_dir, i, build);
|
|
71
|
154 # trackDb.txt track hub entry.
|
|
54
|
155 track_db = c(track_db, paste("track ", hub_name, "_track_", i, sep=""));
|
|
|
156 track_db = c(track_db, "type bigBed");
|
|
|
157 track_db = c(track_db, paste("bigDataUrl", big_data_url, sep=" "));
|
|
77
|
158 track_db = c(track_db, paste("shortLabel", hub_short_label, sep=" "));
|
|
|
159 track_db = c(track_db, paste("longLabel", hub_long_label, sep=" "));
|
|
54
|
160 track_db = c(track_db, paste("priority", i));
|
|
|
161 track_db = c(track_db, "itemRgb on");
|
|
|
162 track_db = c(track_db, "maxItems 100000");
|
|
66
|
163 track_db = c(track_db, paste("color", color, sep=" "));
|
|
54
|
164 track_db = c(track_db, "visibility dense");
|
|
|
165 track_db = c(track_db, "");
|
|
|
166 }
|
|
|
167 return(track_db);
|
|
0
|
168 }
|
|
|
169
|
|
39
|
170 get_big_data_url = function(galaxy_url, encoded_dataset_id, tracks_dir, index, build) {
|
|
56
|
171 track_files <- list.files(path=tracks_dir, pattern="\\.bigbed");
|
|
54
|
172 s <- paste(galaxy_url, 'datasets/', encoded_dataset_id, '/display/myHub/', build, '/', track_files[index], sep="");
|
|
|
173 return(s)
|
|
50
|
174 }
|
|
|
175
|
|
39
|
176 get_track_file_name = function(base_track_file_name, index, ext) {
|
|
54
|
177 track_file_name <- paste(base_track_file_name, index, ext, sep=".");
|
|
|
178 return(track_file_name);
|
|
39
|
179 }
|
|
|
180
|
|
20
|
181 # Create the hub.txt output.
|
|
53
|
182 hub_name_line <- paste("hub ", opt$hub_name, sep="");
|
|
77
|
183 hub_short_label_line <- paste("shortLabel ", opt$hub_short_label, sep="");
|
|
|
184 hub_long_label_line <- paste("longLabel ", opt$hub_long_label, sep="");
|
|
35
|
185 genomes_txt_line <- paste("genomesFile genomes.txt", sep="");
|
|
39
|
186 email_line <- paste("email ", opt$email, sep="");
|
|
77
|
187 contents <- paste(hub_name_line, hub_short_label_line, hub_long_label_line, genomes_txt_line, email_line, sep="\n");
|
|
35
|
188 hub_dir <- paste(opt$output_trackhub_files_path, "/", "myHub", "/", sep="");
|
|
22
|
189 dir.create(hub_dir, showWarnings=FALSE);
|
|
35
|
190 hub_file_path <- paste(hub_dir, "hub.txt", sep="");
|
|
20
|
191 write.table(contents, file=hub_file_path, quote=F, row.names=F, col.names=F);
|
|
|
192
|
|
|
193 # Create the genomes.txt output.
|
|
39
|
194 genome_line <- paste("genome ", opt$build, sep="");
|
|
35
|
195 track_db_line <- paste("trackDb ", opt$build, "/", "trackDb.txt", sep="");
|
|
39
|
196 contents <- paste(genome_line, track_db_line, sep="\n");
|
|
30
|
197 genomes_file_path <- paste(hub_dir, "genomes.txt", sep="");
|
|
20
|
198 write.table(contents, file=genomes_file_path, quote=F, row.names=F, col.names=F);
|
|
|
199
|
|
0
|
200 # Create the tracks.
|
|
68
|
201 heatmap_path <- paste(opt$script_dir, "create_heatmap.R", sep="/");
|
|
|
202 source(heatmap_path);
|
|
36
|
203 tracks_dir <- paste(hub_dir, opt$build, "/", sep="");
|
|
11
|
204 dir.create(tracks_dir, showWarnings=FALSE);
|
|
70
|
205 track_db <- create_track_db(opt$galaxy_url, opt$output_trackhub_id, opt$input_dir_para, opt$input_dir_state, opt$build,
|
|
77
|
206 opt$chrom_len_file, tracks_dir, opt$hub_name, opt$hub_short_label, opt$hub_long_label, opt$state_indexes,
|
|
|
207 opt$state_names, opt$state_colors);
|
|
0
|
208
|
|
20
|
209 # Create the trackDb.txt output.
|
|
35
|
210 track_db_file_path <- paste(tracks_dir, "trackDb.txt", sep="");
|
|
13
|
211 write.table(track_db, file=track_db_file_path, quote=F, row.names=F, col.names=F);
|
|
11
|
212
|
|
20
|
213 # Create the primary HTML dataset.
|
|
78
|
214 create_primary_html(opt$output_trackhub, tracks_dir, opt$build);
|
