changeset 20:696e435f5607 draft

Uploaded
author greg
date Mon, 20 Nov 2017 10:34:37 -0500
parents 8cd8dabdca19
children 5201e615dd7f
files ideas_genome_tracks.R
diffstat 1 files changed, 59 insertions(+), 47 deletions(-) [+]
line wrap: on
line diff
--- a/ideas_genome_tracks.R	Mon Nov 20 10:34:29 2017 -0500
+++ b/ideas_genome_tracks.R	Mon Nov 20 10:34:37 2017 -0500
@@ -6,14 +6,12 @@
 option_list <- list(
     make_option(c("-b", "--build"), action="store", dest="build", help="Genome build"),
     make_option(c("-c", "--chrom_len_file"), action="store", dest="chrom_len_file", help="Chromosome length file"),
-    make_option(c("-d", "--header"),  action="store", dest="header", default=NULL, help="Track header"),
-    make_option(c("-e", "--state_name"),  action="store", dest="state_name", help="State name"),
-    make_option(c("-f", "--hub_id"),  action="store", dest="hub_id", help="Not sure what this is"),
-    make_option(c("-g", "--email"),  action="store", dest="email", help="User email address"),
-    make_option(c("-i", "--cell_info"),  action="store", dest="cell_info", default=NULL, help="Not sure what this is"),
-    make_option(c("-n", "--hub_name"),  action="store", dest="hub_name", default=NULL, help="Not sure what this is"),
+    make_option(c("-e", "--email"),  action="store", dest="email", help="User email address"),
+    make_option(c("-l", "--long_label"), action="store", dest="long_label", help="Hub long label"),
+    make_option(c("-n", "--hub_name"),  action="store", dest="hub_name", default=NULL, help="Hub name without spaces"),
     make_option(c("-p", "--input_dir_para"), action="store", dest="input_dir_para", help="Directory containing .para outputs from IDEAS"),
     make_option(c("-q", "--input_dir_state"), action="store", dest="input_dir_state", help="Directory containing .state outputs from IDEAS"),
+    make_option(c("-s", "--short_label"), action="store", dest="short_label", help="Hub short label"),
     make_option(c("-u", "--output_track_db"),  action="store", dest="output_track_db", help="Output track db file"),
     make_option(c("-w", "--output_trackhub"),  action="store", dest="output_trackhub", help="Output hub file"),
     make_option(c("-x", "--output_trackhub_files_path"),  action="store", dest="output_trackhub_files_path", help="Output hub extra files path")
@@ -42,6 +40,11 @@
     return(state_color);
 }
 
+get_track_file_name = function(base_track_file_name, index, ext) {
+    track_file_name <- paste(base_track_file_name, index, ext, sep=".");
+    return(track_file_name);
+}
+
 get_state_color = function(statemean, markcolor=NULL)
 {   
     if (length(markcolor) == 0) {
@@ -52,10 +55,12 @@
     rg = apply(statemean, 1, range);
     mm = NULL;
     for(i in 1:dim(statemean)[1]) {
-       mm = rbind(mm, (statemean[i,]-rg[1,i])/(rg[2,i]-rg[1,i]+1e-10));
+       mm = rbind(mm, (statemean[i,]-rg[1,i])+1e-10)/(rg[2,i]-rg[1,i]+1e-10));
     }
-    mm = mm^6; 
-    mm = mm / (apply(mm, 1, sum) + 1e-10);
+    mm = mm^6;
+    if(dim(mm)[2] > 1) {
+        mm = mm / (apply(mm, 1, sum) + 1e-10);
+    }
     mycol = mm%*%markcolor;
     s = apply(statemean, 1, max);
     s = (s - min(s)) / (max(s) - min(s) + 1e-10);
@@ -66,7 +71,19 @@
     return(rt);
 }
 
-create_track = function(input_dir_state, chrom_len_file, base_track_file_name, state_color, state_name, header) {
+create_primary_html = function(output_trackhub, tracks_dir) {
+    track_files <- list.files(path=tracks_dir, full.names=TRUE);
+    s <- "<html><head></head><body><ul>"
+        for (i in 1:length(track_files)) {
+        track_file <- track_files[i];
+        s <- paste(s, "<li><a href=\"", track_file, "\">", track_file, "</a></li>");
+    s <- paste(s, "</ul></body></html>"
+    sink(output_trackhub);
+    cat(s);
+    sink();
+}
+
+create_track = function(input_dir_state, chrom_len_file, base_track_file_name, state_color) {
     state_files <- list.files(path=input_dir_state, full.names=TRUE);
     genome_size = read.table(chrom_len_file);
     g = NULL;
@@ -96,14 +113,9 @@
     posst = as.numeric(g1[,3]);
     posed = as.numeric(g1[,4]);
     state = as.matrix(g1[,5:(dim(g1)[2]-1)]);
-    if (length(state_name)==0) {
-        state_name=0:max(state);
-    }
+    state_name = 0:max(state);
     L = dim(g1)[1];
     n = dim(state)[2];
-    if (length(header) > 0) {
-        colnames(g1) = header;
-    }
     cells = as.character(colnames(g1)[5:(dim(g1)[2]-1)]);
     g1 = NULL;
     options(scipen=999);
@@ -118,8 +130,8 @@
         t = c(t, L);
         np = cbind(chr[t], posst[t0], posed[t], tstate[t]);
         x = cbind(np[,1:3], state_name[as.integer(np[,4])+1], 1000, ".", np[,2:3], state_color[as.numeric(np[,4])+1]);
-        track_file_name_bed <- paste(base_track_file_name, i, "bed", sep=".")
-        track_file_name_bigbed <- paste(base_track_file_name, i, "bigbed", sep=".")
+        track_file_name_bed <- get_track_file_name(base_track_file_name, i, "bed");
+        track_file_name_bigbed <- get_track_file_name(base_track_file_name, i, "bigbed");
         write.table(as.matrix(x), track_file_name_bed, quote=F, row.names=F, col.names=F);
         system(paste("bedToBigBed ", track_file_name_bed, chrom_len_file, " ", track_file_name_bigbed));
         system(paste("rm ", track_file_name_bed));
@@ -127,57 +139,57 @@
     return(cells);
 }
 
-create_track_db = function(input_dir_state, chrom_len_file, header, tracks_dir, hub_id, hub_name, state_color, state_name, cell_info) {
-    base_track_file_name <- paste(tracks_dir, hub_id, sep="");
-    cells = create_track(input_dir_state, chrom_len_file, base_track_file_name, state_color, state_name, header);
-    if (length(cell_info) == 0) {
-        cell_info = cbind(cells, cells, cells, "#000000");
-        cell_info = array(cell_info, dim=c(length(cells), 4));
-    }
+create_track_db = function(input_dir_state, chrom_len_file, tracks_dir, hub_name, short_label, long_label, state_color) {
+    base_track_file_name <- paste(tracks_dir, hub_name, sep="");
+    cells = create_track(input_dir_state, chrom_len_file, base_track_file_name, state_color);
+    cell_info = cbind(cells, cells, cells, "#000000");
+    cell_info = array(cell_info, dim=c(length(cells), 4));
     cell_info = as.matrix(cell_info);
     track_db = NULL;
     for(i in 1:length(cells)) {
-        ii = which(cells[i] == cell_info[,1]);
+
+        ii = which(cells[i] == cellinfo[,1]);
         if (length(ii) == 0) {
             next;
         }
         ii = ii[1];
-        track_db = c(track_db, paste("track bigBed", i, sep=""));
-        track_db = c(track_db, paste("priority", ii));
-        track_db = c(track_db, "type bigBed 9 .");
+        # trackDb.txt track entry.
+        track_db = c(track_db, paste(hub_name, "_track_", i, sep=""));
+        track_db = c(track_db, "type bigBed");
+        track_db = c(track_db, paste("bigDataUrl ", get_track_file_name(base_track_file_name, i, "bigbed"));
+        track_db = c(track_db, paste("shortLabel", short_label));
+        track_db = c(track_db, paste("longLabel", long_label);
+        track_db = c(track_db, paste("priority", i));
         track_db = c(track_db, "itemRgb on");
         track_db = c(track_db, "maxItems 100000");
-        track_db = c(track_db, paste("shortLabel", cell_info[ii, 2]));
-        track_db = c(track_db, paste("longLabel", paste(hub_name, cell_info[ii, 3])));
-        track_db = c(track_db, paste("color", paste(c(col2rgb(cell_info[ii, 4])), collapse=",")));
+        track_db = c(track_db, paste("color", paste(c(col2rgb(cellinfo[ii,4])), collapse=",")));
         track_db = c(track_db, "visibility dense");
         track_db = c(track_db, ""); 
     }
     return(track_db);
 }
 
-if (length(opt$hub_name) == 0) {
-    hub_name <- opt$hub_id;
-} else {
-    hub_name <- opt$hub_name;
-}
-
 # Create the color scheme.
 state_color <- create_color_scheme(opt$input_dir_para);
 
+# Create the hub.txt output.
+contents <- c(paste("hub", opt$hub_name), paste("shortLabel",opt$short_label), paste("longLabel", opt$long_label), paste("genomesFile genomes.txt", sep=""), paste("email", opt$email))
+hub_file_path <- paste(opt$output_trackhub_files_path, "/", "hub.txt", sep="")
+write.table(contents, file=hub_file_path, quote=F, row.names=F, col.names=F);
+
+# Create the genomes.txt output.
+contents <- c(paste("genome", opt$build), paste("trackDb ", opt$build, "/", "trackDb.txt", sep=""))
+genomes_file_path <- paste(opt$output_trackhub_files_path, "/", "genomes.txt", sep="")
+write.table(contents, file=genomes_file_path, quote=F, row.names=F, col.names=F);
+
 # Create the tracks.
 tracks_dir <- paste(opt$output_trackhub_files_path, "/", "tracks", "/", sep="");
 dir.create(tracks_dir, showWarnings=FALSE);
-track_db <- create_track_db(opt$input_dir_state, opt$chrom_len_file, opt$header, tracks_dir, opt$hub_id, hub_name, state_color, opt$state_name, opt$cell_info);
+track_db <- create_track_db(opt$input_dir_state, opt$chrom_len_file, tracks_dir, opt$hub_name, opt$short_label, opt$long_label, state_color);
 
-# Write the outputs.
+# Create the trackDb.txt output.
 track_db_file_path <- paste(tracks_dir, "/", paste("trackDb.txt", sep=""), sep="")
 write.table(track_db, file=track_db_file_path, quote=F, row.names=F, col.names=F);
 
-contents <- c(paste("genome", opt$build), paste("trackDb ", opt$build, "/", "trackDb.txt", sep=""))
-genomes_file_path <- paste(opt$output_trackhub_files_path, "/", "genomes.txt", sep="")
-write.table(contents, file=genomes_file_path, quote=F, row.names=F, col.names=F);
-
-contents <- c(paste("hub", opt$hub_id), paste("shortLabel", opt$hub_id), paste("longLabel", hub_name), paste("genomesFile genomes.txt", sep=""), paste("email", opt$email))
-hub_file_path <- paste(opt$output_trackhub_files_path, "/", "hub.txt", sep="")
-write.table(contents, file=hub_file_path, quote=F, row.names=F, col.names=F);
+# Create the primary HTML dataset.
+create_primary_html(opt$output_trackhub, tracks_dir);