Mercurial > repos > iuc > ruvseq
comparison ruvseq.R @ 4:3a083c78896e draft default tip
planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/ruvseq commit 30117fce22f3771c9c0c13e70c3ad14b694de6e2
| author | iuc |
|---|---|
| date | Fri, 21 Apr 2023 14:09:17 +0000 |
| parents | fed9d0350d72 |
| children |
comparison
equal
deleted
inserted
replaced
| 3:d1f7fa5bb3cb | 4:3a083c78896e |
|---|---|
| 1 # setup R error handling to go to stderr | 1 # setup R error handling to go to stderr |
| 2 library("getopt") | 2 library("getopt") |
| 3 options(show.error.messages = F, error = function() { | 3 options(show.error.messages = FALSE, error = function() { |
| 4 cat(geterrmessage(), file = stderr()); q("no", 1, F) | 4 cat(geterrmessage(), file = stderr()) |
| 5 q("no", 1, FALSE) | |
| 5 }) | 6 }) |
| 6 options(stringAsFactors = FALSE, useFancyQuotes = FALSE) | 7 options(stringAsFactors = FALSE, useFancyQuotes = FALSE) |
| 7 | 8 |
| 8 setup_cmdline_options <- function() { | 9 setup_cmdline_options <- function() { |
| 9 args <- commandArgs(trailingOnly = TRUE) | 10 args <- commandArgs(trailingOnly = TRUE) |
| 15 "max_k", "max_k", 1, "double", | 16 "max_k", "max_k", 1, "double", |
| 16 "sample_json", "s", 1, "character", | 17 "sample_json", "s", 1, "character", |
| 17 "plots", "p", 1, "character", | 18 "plots", "p", 1, "character", |
| 18 "header", "H", 0, "logical", | 19 "header", "H", 0, "logical", |
| 19 "txtype", "y", 1, "character", | 20 "txtype", "y", 1, "character", |
| 20 "tx2gene", "x", 1, "character"), # a space-sep tx-to-gene map or GTF file (auto detect .gtf/.GTF) | 21 "tx2gene", "x", 1, "character", # a space-sep tx-to-gene map or GTF file (auto detect .gtf/.GTF) |
| 22 "ruv_ncounts", "ruv_ncounts", 0, "logical"), | |
| 21 byrow = TRUE, ncol = 4) | 23 byrow = TRUE, ncol = 4) |
| 22 | 24 |
| 23 opt <- getopt(spec) | 25 opt <- getopt(spec) |
| 24 # if help was asked for print a friendly message | 26 # if help was asked for print a friendly message |
| 25 # and exit with a non-zero error code | 27 # and exit with a non-zero error code |
| 153 opt <- setup_cmdline_options() | 155 opt <- setup_cmdline_options() |
| 154 alpha <- opt$alpha | 156 alpha <- opt$alpha |
| 155 min_k <- opt$min_k | 157 min_k <- opt$min_k |
| 156 max_k <- opt$max_k | 158 max_k <- opt$max_k |
| 157 min_c <- opt$min_mean_count | 159 min_c <- opt$min_mean_count |
| 160 ruv_ncounts <- ifelse(is.null(opt$ruv_ncounts), FALSE, TRUE) | |
| 158 sample_json <- fromJSON(opt$sample_json) | 161 sample_json <- fromJSON(opt$sample_json) |
| 159 sample_paths <- sample_json$path | 162 sample_paths <- sample_json$path |
| 160 sample_names <- sample_json$label | 163 sample_names <- sample_json$label |
| 161 condition <- as.factor(sample_json$condition) | 164 condition <- as.factor(sample_json$condition) |
| 162 sample_table <- data.frame(samplename = sample_names, filename = sample_paths, condition = condition) | 165 sample_table <- data.frame(samplename = sample_names, filename = sample_paths, condition = condition) |
| 181 set <- result[[name]] | 184 set <- result[[name]] |
| 182 unwanted_variation <- pData(set) | 185 unwanted_variation <- pData(set) |
| 183 df <- data.frame(identifier = rownames(unwanted_variation)) | 186 df <- data.frame(identifier = rownames(unwanted_variation)) |
| 184 df <- cbind(df, unwanted_variation) | 187 df <- cbind(df, unwanted_variation) |
| 185 colnames(df)[2] <- "condition" | 188 colnames(df)[2] <- "condition" |
| 186 write.table(df, file = paste0("batch_effects_", name, ".tabular"), sep = "\t", quote = F, row.names = F) | 189 write.table(df, file = paste0("uv_batch_effects_", name, ".tabular"), sep = "\t", quote = FALSE, row.names = FALSE) |
| 190 if (ruv_ncounts) { | |
| 191 ruvnorm_counts <- normCounts(set) | |
| 192 ruvnorm_df <- data.frame(geneID = rownames(ruvnorm_counts), ruvnorm_counts) | |
| 193 write.table(ruvnorm_df, file = paste0("ruv_norm_counts_", name, ".tabular"), sep = "\t", quote = FALSE, row.names = FALSE) | |
| 194 } | |
| 187 } | 195 } |
| 196 | |
| 188 } | 197 } |
| 189 | 198 |
| 190 # close the plot device | 199 # close the plot device |
| 191 if (!is.null(opt$plots)) { | 200 if (!is.null(opt$plots)) { |
| 192 cat("closing plot device\n") | 201 cat("closing plot device\n") |
