comparison deseq2.R @ 30:8fe98f7094de draft

planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/deseq2 commit 6868b66f73ddbe947986d1a20b546873cbd515a9

29:cd9874cb9019

#   1 "parametric"
#   2 "local"
#   3 "mean"

# setup R error handling to go to stderr
options(show.error.messages = F, error = function() {
  cat(geterrmessage(), file = stderr()); q("no", 1, F)
})

# we need that to not crash galaxy with an UTF8 error on German LC settings.
loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")

  "many_contrasts", "m", 0, "logical",
  "outlier_replace_off", "a", 0, "logical",
  "outlier_filter_off", "b", 0, "logical",
  "auto_mean_filter_off", "c", 0, "logical",
  "beta_prior_off", "d", 0, "logical",
  "alpha_ma", "A", 1, "numeric"
), byrow = TRUE, ncol = 4)
opt <- getopt(spec)

# if help was asked for print a friendly message
# and exit with a non-zero error code

            ## If we have no other information, estimate from raw.
            cat("\nsize factors for samples: no tximport data, no normalization factors, estimating from raw data\n")
            size_factors <- estimateSizeFactorsForMatrix(counts(dds))
        }
    }
    write.table(size_factors, file = opt$sizefactorsfile, sep = "\t", col.names = F, quote = FALSE)
}

apply_batch_factors <- function(dds, batch_factors) {
  rownames(batch_factors) <- batch_factors$identifier
  batch_factors <- subset(batch_factors, select = -c(identifier, condition))

    stop("Batch factor names don't correspond to input sample names, check input files")
  }
  dds_data <- colData(dds)
  # Merge dds_data with batch_factors using indexes, which are sample names
  # Set sort to False, which maintains the order in dds_data
  reordered_batch <- merge(dds_data, batch_factors, by.x = 0, by.y = 0, sort = F)
  batch_factors <- reordered_batch[, ncol(dds_data):ncol(reordered_batch)]
  for (factor in colnames(batch_factors)) {
    dds[[factor]] <- batch_factors[[factor]]
  }
  colnames(dds) <- reordered_batch[, 1]
  return(dds)
}

if (!is.null(opt$batch_factors)) {
  batch_factors <- read.table(opt$batch_factors, sep = "\t", header = T)
  dds <- apply_batch_factors(dds = dds, batch_factors = batch_factors)
  batch_design <- colnames(batch_factors)[-c(1, 2)]
  design_formula <- as.formula(paste("~", paste(c(batch_design, rev(factors)), collapse = " + ")))
  design(dds) <- design_formula
}

  print(sample_table[, -c(1:2), drop = FALSE])
  cat("\ndesign formula:\n")
  print(design_formula)
  cat("\n\n")
  cat(paste(ncol(dds), "samples with counts over", nrow(dds), "genes\n"))
}

# optional outlier behavior
if (is.null(opt$outlier_replace_off)) {
  min_rep <- 7

30:8fe98f7094de

#   1 "parametric"
#   2 "local"
#   3 "mean"

# setup R error handling to go to stderr
options(show.error.messages = FALSE, error = function() {
  cat(geterrmessage(), file = stderr())
  q("no", 1, FALSE)
})

# we need that to not crash galaxy with an UTF8 error on German LC settings.
loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")

  "many_contrasts", "m", 0, "logical",
  "outlier_replace_off", "a", 0, "logical",
  "outlier_filter_off", "b", 0, "logical",
  "auto_mean_filter_off", "c", 0, "logical",
  "beta_prior_off", "d", 0, "logical",
  "alpha_ma", "A", 1, "numeric",
  "prefilter", "P", 0, "logical",
  "prefilter_value", "V", 1, "numeric"
), byrow = TRUE, ncol = 4)
opt <- getopt(spec)

# if help was asked for print a friendly message
# and exit with a non-zero error code

            ## If we have no other information, estimate from raw.
            cat("\nsize factors for samples: no tximport data, no normalization factors, estimating from raw data\n")
            size_factors <- estimateSizeFactorsForMatrix(counts(dds))
        }
    }
    write.table(size_factors, file = opt$sizefactorsfile, sep = "\t", col.names = FALSE, quote = FALSE)
}

apply_batch_factors <- function(dds, batch_factors) {
  rownames(batch_factors) <- batch_factors$identifier
  batch_factors <- subset(batch_factors, select = -c(identifier, condition))

    stop("Batch factor names don't correspond to input sample names, check input files")
  }
  dds_data <- colData(dds)
  # Merge dds_data with batch_factors using indexes, which are sample names
  # Set sort to False, which maintains the order in dds_data
  reordered_batch <- merge(dds_data, batch_factors, by.x = 0, by.y = 0, sort = FALSE)
  batch_factors <- reordered_batch[, ncol(dds_data):ncol(reordered_batch)]
  for (factor in colnames(batch_factors)) {
    dds[[factor]] <- batch_factors[[factor]]
  }
  colnames(dds) <- reordered_batch[, 1]
  return(dds)
}

if (!is.null(opt$batch_factors)) {
  batch_factors <- read.table(opt$batch_factors, sep = "\t", header = TRUE)
  dds <- apply_batch_factors(dds = dds, batch_factors = batch_factors)
  batch_design <- colnames(batch_factors)[-c(1, 2)]
  design_formula <- as.formula(paste("~", paste(c(batch_design, rev(factors)), collapse = " + ")))
  design(dds) <- design_formula
}

  print(sample_table[, -c(1:2), drop = FALSE])
  cat("\ndesign formula:\n")
  print(design_formula)
  cat("\n\n")
  cat(paste(ncol(dds), "samples with counts over", nrow(dds), "genes\n"))
}

# minimal pre-filtering
if (!is.null(opt$prefilter)) {
    keep <- rowSums(counts(dds)) >= opt$prefilter_value
    dds <- dds[keep, ]
}

# optional outlier behavior
if (is.null(opt$outlier_replace_off)) {
  min_rep <- 7