R/write.genes.R

In nickfishbane/DMRtools: Tools for 450K differential methylation analysis

## writes all the gene_assocs found in the dmrcate.output object OR a top-table
## to a simple text list
#' write gene list from results
#'
#' writes all the associated genes found in a \code{dmrcate.output} object or a
#' top-table generated by \code{limma.tt} that satisfy passing a certain threshold
#'
#' if \code{res.obj} is a \code{dmrcate.output} object, it will print all genes with minimum p-value (across probes in the region) less than the supplied threshold.
#' If a threshold is desired, this should be specified in the creation of the
#' object, using \code{pcutoff} in the \code{DMRcate::dmrcate()} function. If
#' \code{res.obj} is a top-table, it will print all gene symbols where the FDR
#' q-value is below the supplied threshold
#'
#' @param res.obj Either a data.frame (limma top-table) with column names "adj.P.Val" and "Gene" or an object of class \code{dmrcate.output} created by the \code{dmrcate} function in the \code{DMRcate} package
#' @param file name (with suggested '.txt' extension) for output of gene list
#' @param thresh Numerical value for p-value cut-off when writing genes from top-table
#'
#' @export
write.genes <- function(res.obj, file, thresh = 0.05){
  if(is(res.obj, 'data.frame')){
    stopifnot(all(c('Gene', 'adj.P.Val') %in% names(res.obj)))
    genes <- subset(res.obj, adj.P.Val < thresh)$Gene
    genes <- unique(unlist(strsplit(genes, ';')))
  } else if(is(res.obj, 'dmrcate.output')){
    cat(sprintf('\n %d DMR regions in the results file \n', nrow(res.obj$results)))
    genes <- subset(res.obj$results, minpval < thresh)$gene_assoc
    genes <- unique(unlist(strsplit(genes, ',')))
    cat(sprintf('%d unique DMR-associated genes identified \n', length(genes)))
  } else stop('res.obj is the wrong kind of format')
  write.table(unique(genes), file = file, col.names = F, row.names = F,
              quote = F, sep = '\n')
}