|
2
|
1 #!/usr/bin/env Rscript
|
|
|
2 argspec <- c("tab.2.cdt.R converts a data matrix to cdt format
|
|
|
3
|
|
|
4 Usage:
|
|
|
5
|
|
|
6 Optional:
|
|
|
7
|
|
|
8 \n\n")
|
|
|
9 args <- commandArgs(TRUE)
|
|
|
10 if ( length( args ) == 1 && args =="--help") {
|
|
|
11 write(argspec, stderr())
|
|
|
12 q();
|
|
|
13 }
|
|
|
14
|
|
|
15 lib.load.quiet <- function( package ) {
|
|
|
16 package <- as.character(substitute(package))
|
|
|
17 suppressPackageStartupMessages( do.call( "library", list( package=package ) ) )
|
|
|
18 }
|
|
|
19
|
|
|
20 lib.load.quiet( getopt )
|
|
|
21 lib.load.quiet( amap )
|
|
|
22 if ( any( c( 'flashClust', 'fastcluster' ) %in% installed.packages() ) ) {
|
|
|
23 if ( 'flashClust' %in% installed.packages() ) {
|
|
|
24 lib.load.quiet( flashClust )
|
|
|
25 } else {
|
|
|
26 if ( 'fastcluster' %in% installed.packages() ) {
|
|
|
27 lib.load.quiet( fastcluster )
|
|
|
28 }
|
|
|
29 }
|
|
|
30 }
|
|
|
31
|
|
|
32 spec <- matrix( c( "dataset", "d", 1, "character",
|
|
|
33 "gen.new.dgram", "g", 2, "character",
|
|
|
34 "output.fname", "o", 2, "character"
|
|
|
35 ),
|
|
|
36 nc=4,
|
|
|
37 byrow=TRUE
|
|
|
38 )
|
|
|
39
|
|
|
40
|
|
|
41 opt <- getopt( spec=spec )
|
|
|
42 if ( is.null( opt$output.report.dir ) ) { opt$output.report.dir <- "report" }
|
|
|
43 if ( is.null( opt$output.fname ) ) { opt$output.fname <- file.path( opt$output.report.dir, paste( "data", opt$output.format, sep="." ) ) }
|
|
|
44 if ( is.null( opt$gen.new.dgram ) ) {
|
|
|
45 opt$gen.new.dgram <- FALSE
|
|
|
46 } else {
|
|
|
47 if ( ! opt$gen.new.dgram %in% c( "no", "yes" ) ) {
|
|
|
48 stop( "invalid input to gen.new.dgram param", opt$gen.new.dgram, "\n" )
|
|
|
49 }
|
|
|
50 ## set to TRUE/FALSE
|
|
|
51 opt$gen.new.dgram <- ( opt$gen.new.dgram == "yes" )
|
|
|
52 }
|
|
|
53
|
|
|
54
|
|
|
55 load( opt$dataset ) ## should load the cl, treecl.res (or partcl.res) and data
|
|
|
56
|
|
|
57 if ( ! exists( 'data' ) ) stop( "No data object in the rdata file provided for", opt$output.format, "format!!\n" )
|
|
|
58 if ( inherits( data, "dist" ) ) stop( "data provided is a distance matrix - not a data matrix. Can't generate TreeView or Tab-delimited files w/distance matrices!\n" )
|
|
|
59
|
|
|
60 ## the rest of this is for the remaining output formats
|
|
|
61 ## pre-set the cluster results for rows & cols to NULL
|
|
|
62 direction <- NULL
|
|
|
63 if ( exists( 'treecl.res' ) ) {
|
|
|
64 cl.res <- treecl.res
|
|
|
65 if ( is.null( treecl.res$dist.method ) ) treecl.res$dist.method <- 'euclidean' # just set it to some stub so that the ctc fn's don't complain
|
|
|
66 } else {
|
|
|
67 if ( exists( 'partcl.res' ) ) {
|
|
|
68 cl.res <- partcl.res
|
|
|
69 }
|
|
|
70 else {
|
|
|
71 stop( 'could not find a valid cluster result to use for primary direction\n' )
|
|
|
72 }
|
|
|
73 }
|
|
|
74
|
|
|
75 if ( all( names( cl ) %in% rownames( data ) ) ) {
|
|
|
76 direction <- "rows"
|
|
|
77 } else if ( all( names( cl ) %in% colnames( data ) ) ) {
|
|
|
78 direction <- "cols"
|
|
|
79 data <- t( data )
|
|
|
80 } else {
|
|
|
81 stop( "Specified cluster result does not come from this data set\n" )
|
|
|
82 }
|
|
|
83
|
|
|
84
|
|
|
85 centroids <- NULL
|
|
|
86 cl <- sort( cl )
|
|
|
87 if ( inherits( cl.res, "kmeans" ) ) {
|
|
|
88 ## already comes pre-calculated for us!!
|
|
|
89 centroids <- cl.res$centers
|
|
|
90 } else {
|
|
|
91 data <- data[ names( cl ), ]
|
|
|
92 cl.list <- unique( cl )
|
|
|
93 cl.list <- lapply( cl.list, function(i) cl[ cl %in% i ] )
|
|
|
94 centroids <- sapply( cl.list,
|
|
|
95 function(x) {
|
|
|
96 return( apply( data[ names(x), ], 2, mean, na.rm=T ) )
|
|
|
97 }
|
|
|
98 )
|
|
|
99 centroids <- t( centroids ) ## get them back to the same number of columns that data has now
|
|
|
100 }
|
|
|
101
|
|
|
102 data <- centroids
|
|
|
103 colnames( data ) <- sapply( 1:max( cl ), function(i) sprintf( "cluster-%02d", i ) )
|
|
|
104
|
|
|
105 if ( opt$gen.new.dgram ) {
|
|
|
106 distance <- 'euclidean'
|
|
|
107 if ( inherits( cl.res, 'hclust' ) ) {
|
|
|
108 distance <- cl.res$dist.method
|
|
|
109 }
|
|
|
110 dmat <- Dist( data, distance )
|
|
|
111 treecl.res <- hclust( dmat )
|
|
|
112 cl <- cutree( treecl.res )
|
|
|
113 }
|
|
|
114
|
|
|
115 if ( direction == "cols" ) {
|
|
|
116 data <- t( data )
|
|
|
117 }
|
|
|
118
|
|
|
119 save( file=opt$output.name, treecl.res, cl, data )
|
