annotate cluster.tools/fix.and.merge.TCGA.sample.IDs.R @ 2:b442996b66ae draft

Uploaded
author peter-waltman
date Wed, 27 Feb 2013 20:17:04 -0500
parents
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
2
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
1 #!/usr/bin/env Rscript
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
2 argspec <- c("fix.and.merge.TCGA.samples.IDs.R takes a clustering from ConsensusClusterPlus and clinical survival data
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
3 and generates a KM-plot, along with the log-rank p-values
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
4
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
5 Usage:
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
6 fix.and.merge.TCGA.samples.IDs.R -d <data.file>
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
7
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
8 \n\n")
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
9 args <- commandArgs(TRUE)
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
10 if ( length( args ) == 1 && args =="--help") {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
11 write(argspec, stderr())
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
12 q();
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
13 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
14
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
15 lib.load.quiet <- function( package ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
16 package <- as.character(substitute(package))
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
17 suppressPackageStartupMessages( do.call( "library", list( package=package ) ) )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
18 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
19 lib.load.quiet(getopt)
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
20
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
21 spec <- matrix( c( "data.fname", "d", 1, "character",
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
22 "num.components", "n", 2, "integer",
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
23 "remove.normals", "r", 0, "logical",
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
24 "output.fname", "o", 2, "character"
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
25 ),
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
26 nc=4,
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
27 byrow=TRUE
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
28 )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
29
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
30 opt <- getopt( spec=spec )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
31
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
32 data <- as.matrix( read.delim( opt$data.fname, row.names=1, check.names=FALSE ) )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
33 if ( is.null( opt$num.components ) ) { opt$num.components <- 3 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
34 if ( is.null( opt$remove.normals ) ) { opt$remove.normals <- FALSE }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
35 if ( is.null( opt$output.fname ) ) { opt$output.fname <- paste( "sample.IDs.updated", basename( opt$data.fname ), sep="." ) }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
36
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
37 if ( opt$num.components < 3 ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
38 err.msg <- "Minimum number of barcode components that can be used is 3\n"
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
39 cat( err.msg, file=opt$output.fname )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
40 stop( err.msg )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
41 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
42
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
43 remove.periods.from.ids <- function( ids ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
44 return( gsub( "\\.", "-", ids ) )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
45 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
46
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
47
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
48 reformat.ids <- function( ids,
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
49 num.components=3 ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
50 return( sapply( strsplit( ids, "-" ), function(x) paste( x[1:num.components], collapse="-" ) ) )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
51 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
52
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
53
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
54 merge.cols <- function( mat,
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
55 samp.ids ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
56
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
57 if ( ! any( duplicated( samp.ids ) ) ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
58 colnames( mat ) <- samp.ids
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
59 return( mat )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
60 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
61
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
62 dupes <- unique( samp.ids[ duplicated( samp.ids ) ] )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
63 uniqs <- samp.ids[ ! samp.ids %in% dupes ]
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
64
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
65 uniq.mat <- mat[ , ( samp.ids %in% uniqs ), drop=FALSE ]
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
66 colnames( uniq.mat ) <- uniqs
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
67
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
68 for ( dup in dupes ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
69 dup.mat <- apply( mat[, ( samp.ids %in% dup ), drop=FALSE],
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
70 1,
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
71 mean,
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
72 na.rm=TRUE )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
73
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
74 uniq.mat <- cbind( uniq.mat, dup.mat )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
75 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
76 colnames( uniq.mat ) <- c( uniqs, dupes )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
77 return( uniq.mat )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
78 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
79
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
80
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
81 cnames <- colnames( data )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
82 rnames <- rownames( data )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
83
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
84 transpose.back <- FALSE
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
85
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
86 if ( all( grepl( "^TCGA", rnames ) ) ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
87 data <- t( data )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
88 transpose.back <- TRUE
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
89 } else {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
90 if ( ! all( grepl( "^TCGA", cnames ) ) ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
91 err.msg <- "can't find any TCGA samples listed in this matrix. If columns are samples, all columns must be a TCGA sample ID. Same if rows are samples.\n"
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
92 cat( err.msg, file=opt$output.fname )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
93 stop( err.msg )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
94 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
95 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
96
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
97 cnames <- remove.periods.from.ids( colnames( data ) )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
98 nelts <- as.numeric( names( table( as.factor( sapply( strsplit( cnames, "-" ), function(x) length(x ) ) ) ) ) )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
99 if ( length( nelts ) > 1 ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
100 err.msg <- "Error: Inconsistent TCGA sample barcodes used. Have found ID with different numbers of components in the barcodes used\n"
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
101 cat( err.msg, file=opt$output.fname )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
102 stop( err.msg )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
103 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
104
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
105 if ( opt$remove.normals ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
106 if ( nelts > 3 ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
107 normals <- grepl( "^TCGA-..-....-1", cnames )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
108 data <- data[ , (! normals ), drop=FALSE ]
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
109 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
110 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
111
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
112 if ( opt$num.components < nelts ) {
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
113 cnames <- reformat.ids( ids=cnames, num.components=opt$num.components )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
114 data <- merge.cols( data, cnames )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
115 }
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
116
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
117 if ( transpose.back ) data <- t( data )
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
118
b442996b66ae Uploaded
peter-waltman
parents:
diff changeset
119 write.table( data, opt$output.fname, sep="\t", quote=FALSE, col.names=NA )