7
|
1 #!/usr/bin/env Rscript
|
|
2 argspec <- c("tab.2.cdt.R converts a data matrix to cdt format
|
|
3
|
|
4 Usage:
|
|
5 tab.2.cdt.R -d <data.file>
|
|
6 Optional:
|
|
7 -o <output_file>
|
|
8 \n\n")
|
|
9 args <- commandArgs(TRUE)
|
|
10 if ( length( args ) == 1 && args =="--help") {
|
|
11 write(argspec, stderr())
|
|
12 q();
|
|
13 }
|
|
14
|
|
15 ## some helper fn's
|
|
16 write.2.tab <- function( mat,
|
|
17 fname ) {
|
|
18 mat <- rbind( colnames( mat ), mat )
|
|
19 mat <- cbind( c( "ID", rownames( mat )[-1] ),
|
|
20 mat )
|
|
21 write.table( mat, fname, sep="\t", row.names=FALSE, col.names=FALSE, quote=FALSE )
|
|
22 }
|
|
23
|
|
24 lib.load.quiet <- function( package ) {
|
|
25 package <- as.character(substitute(package))
|
|
26 suppressPackageStartupMessages( do.call( "library", list( package=package ) ) )
|
|
27 }
|
|
28
|
|
29 lib.load.quiet( getopt )
|
|
30 lib.load.quiet( ctc )
|
|
31 if ( any( c( 'flashClust', 'fastcluster' ) %in% installed.packages() ) ) {
|
|
32 if ( 'flashClust' %in% installed.packages() ) {
|
|
33 lib.load.quiet( flashClust )
|
|
34 } else {
|
|
35 if ( 'fastcluster' %in% installed.packages() ) {
|
|
36 lib.load.quiet( fastcluster )
|
|
37 }
|
|
38 }
|
|
39 }
|
|
40
|
|
41 spec <- matrix( c( "dataset", "d", 1, "character",
|
|
42 "num.k", "k", 2, "character",
|
|
43 "output.fname", "o", 2, "character"
|
|
44 ),
|
|
45 nc=4,
|
|
46 byrow=TRUE
|
|
47 )
|
|
48
|
|
49
|
|
50 opt <- getopt( spec=spec )
|
|
51 if ( is.null( opt$output.fname ) ) { opt$output.fname <- file.path( opt$output.report.dir, paste( "data", opt$output.format, sep="." ) ) }
|
|
52 if ( is.null( opt$num.k ) ) {
|
|
53 opt$num.k <- -1
|
|
54 } else {
|
|
55 num.k <- as.integer( eval( parse( text=paste( "c(",gsub( "-", ":", gsub( ", |,", ",", opt$num.k ) ), ")" ) ) ) )
|
|
56 num.k <- num.k[ ! is.na( num.k ) ]
|
|
57 if ( length( opt$num.k ) == 0 ) stop( 'invalid input for k_range specified:', opt$num.k, "\n" )
|
|
58
|
|
59 num.k <- num.k[ ! num.k %in% 1 ] # strip out a k==1 since that doesn't make any sense
|
|
60 opt$num.k <- num.k; rm( num.k )
|
|
61 }
|
|
62
|
|
63 load( opt$dataset )
|
|
64 ## if this is a clustering result w/cluster assignments ('raw' CCPLUS does not)
|
|
65 if ( exists( 'cl' ) ) {
|
|
66 k <- max( as.numeric( cl ) )
|
|
67 cl <- matrix( cl, nc=1, dimnames=list( names(cl), k ) )
|
|
68 if ( (length(opt$num.k)==1) && (opt$num.k == -1 ) ) opt$num.k <- k
|
|
69
|
|
70 ## if this is a one-off to produce a phenotype for the number of clusters that the user originally proposed
|
|
71 if ( !opt$num.k[1] %in% c( -1, k ) ) {
|
|
72
|
|
73 if ( exists( 'partcl.res' ) || exists( 'select.result' ) ) {
|
|
74
|
|
75 if ( exists( 'partcl.res' ) ) {
|
|
76 warning( 'The k_range value(s) specified are:',
|
|
77 opt$num.k,
|
|
78 "however k_range vals can not specify alternate k values for partition clusters. Using the K value that corresponds to this result instead\n" )
|
|
79 } else {
|
|
80 warning( 'The k_range value(s) specified are:',
|
|
81 opt$num.k,
|
|
82 "however k_range vals can not specify alternate k values for specific cluster results from CCPLUS (i.e. those from the Select K or Extract tools). To get alternate K values, re-run the dichotomizer on the 'raw' CCPLUS results. Using the K value that corresponds to this result instead\n" )
|
|
83 }
|
|
84
|
|
85 opt$num.k <- k
|
|
86 cl <- matrix( cl, nr=1, dimnames=list( k, names(cl) ) )
|
|
87 } else {
|
|
88 ## handle if this is a hclust result
|
|
89
|
|
90 opt$num.k <- opt$num.k[ opt$num.k < length( cl ) ]
|
|
91 cl.samps <- rownames( cl )
|
|
92 cl <- sapply( opt$num.k, function(i) cutree( treecl.res, i )[ cl.samps ] )
|
|
93 colnames( cl ) <- opt$num.k
|
|
94 }
|
|
95 }
|
|
96 } else if ( exists( 'results' ) ) {
|
|
97 ## handle if this is a ccplus-raw result
|
|
98 opt$num.k <- opt$num.k[ opt$num.k <= length( results ) ]
|
|
99 cl <- sapply( results[ opt$num.k ], '[[', 'consensusClass' )
|
|
100 colnames( cl ) <- opt$num.k
|
|
101 }
|
|
102
|
|
103 pheno.mat <- lapply( 1:ncol(cl),
|
|
104 function(i) {
|
|
105 x <- cl[,i]
|
|
106 cls <- ks <- sort( unique(x) )
|
|
107 cls <- sapply( cls, function(y) as.numeric( x %in% y ) )
|
|
108 colnames(cls) <- paste( "CLeq", ks, sep="" )
|
|
109 rownames(cls) <- names(x)
|
|
110 return(cls)
|
|
111 }
|
|
112 )
|
|
113 names( pheno.mat ) <- opt$num.k
|
|
114
|
|
115 final.mat <- matrix( NA, nc=0, nrow=nrow(cl), dimnames=list( names(cl), NULL ) )
|
|
116 for ( i in names( pheno.mat ) ) {
|
|
117 colnames( pheno.mat[[i]] ) <- paste( "Keq", i, "_", colnames( pheno.mat[[i]] ), sep="" )
|
|
118 final.mat <- cbind( final.mat, pheno.mat[[i]] )
|
|
119 }
|
|
120
|
|
121 write.2.tab( final.mat, opt$output.fname )
|