R/Step4.R

#' Run limma on output of Step3
#' @param REVIEW3 output from step3
#' @param CONTRASTS output from step3
#' @param datadir path to GEOtemp, e.g.: "GEOtemp/" backslash is required
#' @param resultsdir path to output directory, e.g.: "results/" backslash is required
limma_fit <- function(REVIEW3, CONTRASTS, datadir, resultsdir){
  ## Fit Design (Fit 1)

  DT <- fread(REVIEW3, colClasses = c("RefDye" = "character", "Comment" = "character", "Group" = "character", "Batch" = "character",
                                      "Node" = "character", "NodeFunction" = "character"))
  DT <- DT[GSE != "", 1:17]
  setkey(DT, GSE)
  contrasts <- fread(CONTRASTS)
  contrasts <- split(contrasts, by = c("GSE", "BioSampName"))

  fit1 <- lapply(contrasts, try2getFit1)


  ## Fit Contrasts (Fit 2)

  fit2 <- lapply(fit1, function(x) try2getFit2(x$Fit, x$Contrast))
  ok <- sapply(fit2, allOK)
  results <- fit2[ok]
  results <- results[!sapply(results, is.null)]

  for(gse in names(results)) {
    write.fit(results[[gse]], results = NULL, file = paste0(resultsdir, gse, "_fit.txt"), digits = 10,
              adjust = "fdr", method = "separate", F.adjust = "none", sep="\t")
  }
}
axelmuller/geometric2 documentation built on May 25, 2019, 6:26 p.m.