# HG changeset patch # User greg # Date 1513363292 18000 # Node ID d7ab58fe039f5abf30be375a4db251f214ee1fa5 # Parent 8b1c49f2ecaaa6081756bd0e521e372c6cad58c9 Uploaded diff -r 8b1c49f2ecaa -r d7ab58fe039f ideas_genome_tracks.R --- a/ideas_genome_tracks.R Thu Dec 14 14:49:54 2017 -0500 +++ b/ideas_genome_tracks.R Fri Dec 15 13:41:32 2017 -0500 @@ -5,17 +5,17 @@ suppressPackageStartupMessages(library("viridisLite")) option_list <- list( - make_option(c("--build"), action="store", dest="build", help="Genome build"), - make_option(c("--chrom_len_file"), action="store", dest="chrom_len_file", help="Chromosome length file"), - make_option(c("--email"), action="store", dest="email", help="User email address"), - make_option(c("--galaxy_url"), action="store", dest="galaxy_url", help="Galaxy instance base URL"), - make_option(c("--hub_name"), action="store", dest="hub_name", default=NULL, help="Hub name without spaces"), - make_option(c("--input_dir_state"), action="store", dest="input_dir_state", help="Directory containing .state outputs from IDEAS"), - make_option(c("--long_label"), action="store", dest="long_label", help="Hub long label"), - make_option(c("--output_trackhub"), action="store", dest="output_trackhub", help="Output hub file"), - make_option(c("--output_trackhub_files_path"), action="store", dest="output_trackhub_files_path", help="Output hub extra files path"), - make_option(c("--output_trackhub_id"), action="store", dest="output_trackhub_id", help="Encoded output_trackhub dataset id"), - make_option(c("--short_label"), action="store", dest="short_label", help="Hub short label") + make_option(c("--build"), action="store", dest="build", help="Genome build"), + make_option(c("--chrom_len_file"), action="store", dest="chrom_len_file", help="Chromosome length file"), + make_option(c("--email"), action="store", dest="email", help="User email address"), + make_option(c("--galaxy_url"), action="store", dest="galaxy_url", help="Galaxy instance base URL"), + make_option(c("--hub_name"), action="store", dest="hub_name", default=NULL, help="Hub name without spaces"), + make_option(c("--input_dir_state"), action="store", dest="input_dir_state", help="Directory containing .state outputs from IDEAS"), + make_option(c("--long_label"), action="store", dest="long_label", help="Hub long label"), + make_option(c("--output_trackhub"), action="store", dest="output_trackhub", help="Output hub file"), + make_option(c("--output_trackhub_files_path"), action="store", dest="output_trackhub_files_path", help="Output hub extra files path"), + make_option(c("--output_trackhub_id"), action="store", dest="output_trackhub_id", help="Encoded output_trackhub dataset id"), + make_option(c("--short_label"), action="store", dest="short_label", help="Hub short label") ) parser <- OptionParser(usage="%prog [options] file", option_list=option_list) @@ -23,120 +23,120 @@ opt <- args$options create_primary_html = function(output_trackhub, tracks_dir, build) { - track_files <- list.files(path=tracks_dir); - s <- paste('', sep="\n"); - s <- paste(s, '

Contents of directory ~/myHub/', build, ' required by UCSC TrackHub

\n', sep=""); - s <- paste(s, '\n\n', sep=""); - cat(s, file=output_trackhub); + track_files <- list.files(path=tracks_dir); + s <- paste('', sep="\n"); + s <- paste(s, '

Contents of directory ~/myHub/', build, ' required by UCSC TrackHub

\n', sep=""); + s <- paste(s, '\n\n', sep=""); + cat(s, file=output_trackhub); } create_track = function(input_dir_state, chrom_len_file, base_track_file_name) { - # Create everythin needed, including the bigbed file, - # to render the tracks within the UCSC track hub. - state_files <- list.files(path=input_dir_state, full.names=TRUE); - genome_size = read.table(chrom_len_file); - g = NULL; - for(i in state_files) { - tg = as.matrix(fread(i)); - t = NULL; - for(j in 1:dim(genome_size)[1]) { - t = c(t, which(tg[,2]==as.character(genome_size[j, 1]) & as.numeric(tg[,4]) > as.numeric(genome_size[j, 2]))); - } - if (length(t) > 0) { - tg = tg[-t,]; - } - t = which(is.na(match(tg[,2], genome_size[,1]))==T); - if (length(t)>0) { - tg = tg[-t,]; - } - g = rbind(g, tg); - } - uchr = sort(unique(as.character(g[,2]))); - g1 = NULL; - for(i in uchr) { - t = which(g[,2]==i); - g1 = rbind(g1, g[t[order(as.integer(g[t, 3]))],]); - } - g = NULL; - chr = as.character(g1[,2]); - posst = as.numeric(g1[,3]); - posed = as.numeric(g1[,4]); - state = as.matrix(g1[,5:(dim(g1)[2]-1)]); - state_name = 0:max(state); - L = dim(g1)[1]; - n = dim(state)[2]; - cells = as.character(colnames(g1)[5:(dim(g1)[2]-1)]); - options(scipen=999); - tt = which(chr[2:L]!=chr[2:L-1]); - tt = c(tt, which(posst[2:L]!=posed[2:L-1])); - tt = sort(unique(tt)); - for(i in 1:n) { - tstate = state[,i]; - t = c(tt, which(tstate[2:L]!=tstate[2:L-1])); - t = sort(unique(t)); - t0 = c(0, t) + 1; - t = c(t, L); - np = cbind(chr[t], posst[t0], posed[t], tstate[t]); - 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"); - x = cbind(np[, 1:3], state_name[as.integer(np[,4])+1], 1000, ".", np[,2:3]); - 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)); - } - return(cells); + # Create everythin needed, including the bigbed file, + # to render the tracks within the UCSC track hub. + state_files <- list.files(path=input_dir_state, full.names=TRUE); + genome_size = read.table(chrom_len_file); + g = NULL; + for(i in state_files) { + tg = as.matrix(fread(i)); + t = NULL; + for(j in 1:dim(genome_size)[1]) { + t = c(t, which(tg[,2]==as.character(genome_size[j, 1]) & as.numeric(tg[,4]) > as.numeric(genome_size[j, 2]))); + } + if (length(t) > 0) { + tg = tg[-t,]; + } + t = which(is.na(match(tg[,2], genome_size[,1]))==T); + if (length(t)>0) { + tg = tg[-t,]; + } + g = rbind(g, tg); + } + uchr = sort(unique(as.character(g[,2]))); + g1 = NULL; + for(i in uchr) { + t = which(g[,2]==i); + g1 = rbind(g1, g[t[order(as.integer(g[t, 3]))],]); + } + g = NULL; + chr = as.character(g1[,2]); + posst = as.numeric(g1[,3]); + posed = as.numeric(g1[,4]); + state = as.matrix(g1[,5:(dim(g1)[2]-1)]); + state_name = 0:max(state); + L = dim(g1)[1]; + n = dim(state)[2]; + cells = as.character(colnames(g1)[5:(dim(g1)[2]-1)]); + options(scipen=999); + tt = which(chr[2:L]!=chr[2:L-1]); + tt = c(tt, which(posst[2:L]!=posed[2:L-1])); + tt = sort(unique(tt)); + for(i in 1:n) { + tstate = state[,i]; + t = c(tt, which(tstate[2:L]!=tstate[2:L-1])); + t = sort(unique(t)); + t0 = c(0, t) + 1; + t = c(t, L); + np = cbind(chr[t], posst[t0], posed[t], tstate[t]); + 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, "bb"); + x = cbind(np[, 1:3], state_name[as.integer(np[,4])+1], 1000, ".", np[,2:3]); + 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)); + } + return(cells); } create_track_db = function(galaxy_url, encoded_dataset_id, input_dir_state, build, chrom_len_file, tracks_dir, hub_name, short_label, long_label) { - # Create a trackDb.txt file that includes each state. - base_track_file_name <- paste(tracks_dir, hub_name, sep=""); - cells = create_track(input_dir_state, chrom_len_file, base_track_file_name); - # Create a a character vector of 1024 viridis color hex codes. - viridis_vector <- viridis(1024, alpha=1, begin=0, end=1, direction=1, option="D"); - colors_used <- vector(); - track_db = NULL; - for (i in 1:length(cells)) { - # Generate a random number between 1 and 1024 as - # the viridis_vector index for the next state color. - color_index <- sample(1:1024, 1) - # Make sure the color was not previously chosen. - while(is.element(color_index, colors_used)) { - color_index <- sample(1:1024, 1) - } - # Append the color to our list of chosen colors. - append(colors_used, color_index) - # Get the hex code from viridis_vector. - viridis_color_hex_code <- viridis_vector[color_index] - big_data_url <- get_big_data_url(galaxy_url, encoded_dataset_id, tracks_dir, i, build); - # trackDb.txt track entry. - track_db = c(track_db, paste("track ", hub_name, "_track_", i, sep="")); - track_db = c(track_db, "type bigBed"); - track_db = c(track_db, paste("bigDataUrl", big_data_url, sep=" ")); - track_db = c(track_db, paste("shortLabel", short_label, sep=" ")); - track_db = c(track_db, paste("longLabel", long_label, sep=" ")); - 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("color", paste(c(col2rgb(viridis_color_hex_code)), collapse=","), sep=" ")); - track_db = c(track_db, "visibility dense"); - track_db = c(track_db, ""); - } - return(track_db); + # Create a trackDb.txt file that includes each state. + base_track_file_name <- paste(tracks_dir, hub_name, sep=""); + cells = create_track(input_dir_state, chrom_len_file, base_track_file_name); + # Create a a character vector of 1024 viridis color hex codes. + viridis_vector <- viridis(1024, alpha=1, begin=0, end=1, direction=1, option="D"); + colors_used <- vector(); + track_db = NULL; + for (i in 1:length(cells)) { + # Generate a random number between 1 and 1024 as + # the viridis_vector index for the next state color. + color_index <- sample(1:1024, 1) + # Make sure the color was not previously chosen. + while(is.element(color_index, colors_used)) { + color_index <- sample(1:1024, 1) + } + # Append the color to our list of chosen colors. + append(colors_used, color_index) + # Get the hex code from viridis_vector. + viridis_color_hex_code <- viridis_vector[color_index] + big_data_url <- get_big_data_url(galaxy_url, encoded_dataset_id, tracks_dir, i, build); + # trackDb.txt track entry. + track_db = c(track_db, paste("track ", hub_name, "_track_", i, sep="")); + track_db = c(track_db, "type bigBed"); + track_db = c(track_db, paste("bigDataUrl", big_data_url, sep=" ")); + track_db = c(track_db, paste("shortLabel", short_label, sep=" ")); + track_db = c(track_db, paste("longLabel", long_label, sep=" ")); + 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("color", paste(c(col2rgb(viridis_color_hex_code)), collapse=","), sep=" ")); + track_db = c(track_db, "visibility dense"); + track_db = c(track_db, ""); + } + return(track_db); } get_big_data_url = function(galaxy_url, encoded_dataset_id, tracks_dir, index, build) { - track_files <- list.files(path=tracks_dir, pattern="\\.bigbed"); - s <- paste(galaxy_url, 'datasets/', encoded_dataset_id, '/display/myHub/', build, '/', track_files[index], sep=""); - return(s) + track_files <- list.files(path=tracks_dir, pattern="\\.bigbed"); + s <- paste(galaxy_url, 'datasets/', encoded_dataset_id, '/display/myHub/', build, '/', track_files[index], sep=""); + return(s) } 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); + track_file_name <- paste(base_track_file_name, index, ext, sep="."); + return(track_file_name); } # Create the hub.txt output. @@ -168,4 +168,4 @@ write.table(track_db, file=track_db_file_path, quote=F, row.names=F, col.names=F); # Create the primary HTML dataset. -create_primary_html(opt$output_trackhub, tracks_dir, opt$build); \ No newline at end of file +create_primary_html(opt$output_trackhub, tracks_dir, opt$build);