Commit 76c5dd98 authored by Holger Brandl's avatar Holger Brandl

require contrast attribute to come first in design

parent 5b5cddd2
......@@ -45,8 +45,7 @@ gene_info_file = opts$gene_info
assert(is.null(gene_info_file) || file.exists(gene_info_file), "invalid gene_info_file")
designFormula = opts$design
assert(str_detect(designFormula, ".*condition$"))
assert(str_detect(designFormula, "^condition.*")) ## make sure that the condition comes before all batch factors
results_prefix = if (str_length(opts$out) > 0) opts$out else "" # used by add_prefix
......@@ -194,11 +193,15 @@ orderMatcheExpDesign = data_frame(replicate = colnames(expMatrix)) %>%
right_join(expDesign, by = "replicate") %>%
arrange(col_index)
## build design matrix
#A key strength of limma’s linear modelling approach, is the ability accommodate arbitrary experimental complexity. Simple designs, such as the one in this workflow, with cell type and batch, through to more complicated factorial designs and models with interaction terms can be handled relatively easily
#' Build design matrix
#' > A key strength of limma’s linear modelling approach, is the ability accommodate arbitrary experimental complexity. Simple designs, such as the one in this workflow, with cell type and batch, through to more complicated factorial designs and models with interaction terms can be handled relatively easily
#'
#' Make sure that non of the batch-factors is confounded with treatment (condition). See https://support.bioconductor.org/p/39385/ for a discussion
#' References
#' * https://f1000research.com/articles/5-1408/v1
# design <- orderMatcheExpDesign %$% model.matrix(~ 0 + condition)
design <- orderMatcheExpDesign %$% model.matrix(formula(as.formula(paste("~0+", designFormula))))
# design <- orderMatcheExpDesign %$% model.matrix(~ 0 + condition + prep_day)
design = orderMatcheExpDesign %$% model.matrix(formula(as.formula(paste("~0+", designFormula))))
rownames(design) <- orderMatcheExpDesign$replicate
#design <- model.matrix(~0+group+lane)
......@@ -333,7 +336,8 @@ deResults %<>% left_join(sampleMeans)
#' apply the hit criterion
deResults %>% ggplot(aes(adj_p_val)) + geom_histogram()
deResults %>% ggplot(aes(p_value)) + geom_histogram(binwidth=.01)
deResults %>% ggplot(aes(adj_p_val)) + geom_histogram(binwidth=.01)
# report hit criterion
#+ results='asis'
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment