view cluster.tools/extract.cons.cluster.from.result.R @ 2:b442996b66ae draft

Uploaded

#!/usr/bin/env Rscript
# Extract k cluster assignment from consensus clustering result Script by Peter Waltman
# Nov. 12, 2012
# License under Creative Commons Attribution 3.0 Unported (CC BY 3.0)
#
##usage, options and doc goes here
argspec <- c("galaxy.extract.cons.clustering.from.result.R extracts a cluster assignment for a specified
value of K from a specified consensus cluster result file.


        Usage: 
                galaxy.extract.cons.cluster.from.result.R -r <results_file> -k <k_select>
        Optional:
                -o consensus class output file # tab-delimitted file format
\n\n")
args <- commandArgs(TRUE)
if ( length( args ) == 1 && args =="--help") { 
  write( argspec, stderr() )
  q();
}

library(getopt)

spec <- matrix( c( "results.file",   "r", 1, "character",
                   "k.select",       "k", 1, "integer",
                   "cluster.class.out",   "o", 2, "character",
                   "output.select.rdata", "d", 2, "character"
                  ),
                nc=4,
                byrow=T
               )
opt <- getopt( spec=spec )
if ( is.null( opt$output.select.rdata ) ) { opt$output.select.rdata <- "select.RData" }
##if ( is.null( opt$cluster.class.out) ) { opt$cluster.class.out <- "select.cls" }

load( opt$results.file )
cons.matrices <- lapply( results[ 2:length(results) ], '[[', 'consensusMatrix' )
cls <- lapply( results[ 2:length(results) ], '[[', 'consensusClass' )
names( cons.matrices ) <- names( cls ) <- 2:length( results )

ch.k.select <- as.character( opt$k.select )
if ( ch.k.select %in% names( cls ) ) {
  ## get the consensusClass file that's associated with the k.select
  cl <- cls[[ ch.k.select  ]] 

  if ( ! is.null( opt$cluster.class.out ) ) {
    cl <- cbind( names( cl ), as.integer(cl) )
    colnames( cl ) <- c( "ID", "class" )
    write.table( cl, opt$cluster.class.out, sep="\t", row.names=FALSE, quote=FALSE )
  } else if ( ! is.null( opt$output.select.rdata ) ) {
    ## re-order the samples to follow the cluster assignment
    treecl.res <- results[[ opt$k.select ]]$consensusTree
    select.result <- results[[ opt$k.select ]]
    save( file=opt$output.select.rdata, treecl.res, cl, select.result, data )
  } else {
    stop( 'no valid output format specified\n' )
  }
} else {
  out.string <- paste( "choice of k =", ch.k.select, "not available in this result file. Max k = ", max( as.numeric( names(cls) ) ), "\n" )
  cat( out.string, file=opt$cluster.class.out )
}