|
88
|
1 #!/usr/bin/Rscript
|
|
|
2
|
|
|
3 ## --------------------
|
|
|
4 ## prints all arguments as msg
|
|
|
5 ## --------------------
|
|
|
6 catMsg <- function( msg=c() ){
|
|
|
7 cat( MAIN_NAME, paste( msg, collapse="" ), "\n", sep='')
|
|
|
8 }
|
|
|
9
|
|
|
10 ## ==================================================
|
|
|
11 ## Start of analysis
|
|
|
12 ## ==================================================
|
|
|
13 MAIN_NAME <- '[INFO] '
|
|
|
14 catMsg( "Starting QDNAseq-export wrapper" )
|
|
|
15
|
|
|
16 ## supress msg to allow R to finish with non-error msg
|
|
|
17 catMsg( "Loading R libraries" )
|
|
|
18 suppressWarnings( suppressMessages( library( QDNAseq, quietly = TRUE ) ) )
|
|
|
19 suppressWarnings( suppressMessages( library( CGHcall, quietly = TRUE ) ) )
|
|
|
20 suppressWarnings( suppressMessages( library( CGHregions, quietly = TRUE ) ) )
|
|
|
21
|
|
|
22 ## only one param: the tmp config file
|
|
|
23 cmdLineArgs <- commandArgs(TRUE)
|
|
|
24 config <- cmdLineArgs[1]
|
|
|
25
|
|
|
26 ## sourcing the config file will load all input params
|
|
|
27 ## many variables are imported via sourced "config"
|
|
|
28 source( config ) # outputFile, outputName
|
|
|
29
|
|
|
30 systemUser <- system("whoami",T)
|
|
|
31 qdnaseqVersion <- packageDescription( "QDNAseq" )$Version
|
|
|
32 rVersion <- R.version.string
|
|
|
33 rPath <- R.home()
|
|
|
34 catMsg( c("QDNAseq version ", qdnaseqVersion) )
|
|
|
35 catMsg( c( rVersion ) )
|
|
|
36
|
|
|
37 qdnaseqObject <- readRDS( rdsFilePath )
|
|
|
38 sampleNames <- sampleNames( qdnaseqObject )
|
|
|
39
|
|
|
40 ## sanity checks
|
|
|
41 elements <- assayDataElementNames(qdnaseqObject)
|
|
|
42 if ( ! "calls" %in% elements ) stop( "No calls present in object, regioning with CGHregions only work on called data" )
|
|
|
43 if ( length(sampleNames) < 2 ) stop( "Object contains too few samples, regioning with CGHregions only works with at least two samples" )
|
|
|
44
|
|
|
45 ## analysis
|
|
|
46 cgh <- makeCgh( qdnaseqObject, filter=TRUE )
|
|
|
47 regions <- CGHregions( cgh, averror=0.00001 )
|
|
|
48 outputData <- cbind( regions@featureData@data, regions@assayData$regions )
|
|
|
49
|
|
|
50 ## output
|
|
|
51 write.table( outputData, outputFile, sep="\t", quote=FALSE, row.names=FALSE )
|
|
|
52
|
|
|
53 ## tell galaxy all seems ok
|
|
|
54 q(status=0)
|
