diff roi_metrics.Rscript @ 0:61d9bdb6d519 draft

Uploaded

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/roi_metrics.Rscript	Thu Apr 18 08:03:38 2013 -0400
@@ -0,0 +1,103 @@
+#!/usr/bin/env Rscript
+
+# Compute metrics (all except PCA) for ROI files.
+#
+# USAGE: roi_metrics.Rscript -i IN.roi -o OUT.roi
+
+# TODO(holtgrew): Bernd, can you give the English names of the metrics below in the comments?
+
+# Author: Bernd Jagla <[e-mail]>
+# Author: Manuel Holtgrewe <[e-mail]>
+
+library('getopt')  # parsing args
+library('ngsroi')  # ROI I/O
+library('boot')    # k3.linear()
+
+# ----------------------------------------------------------------------------
+# Command Line Parsing
+# ----------------------------------------------------------------------------
+
+# Parse command line options.
+spec = matrix(c(
+  "help",     "h", 0, "logical",
+  "infile",   "i", 1, "character",
+  "outfile",  "o", 1, "character"
+), byrow=TRUE, ncol=4);
+opt = getopt(spec);
+
+# If help was asked for print a friendly message and exit with a non-zero
+# error code.
+if (!is.null(opt$help) )
+{
+  cat(getopt(spec, usage=TRUE));
+  q(status=1);
+}
+
+# Check arguments.
+if (is.null(opt$infile))
+{
+    write("-i/--infile missing.", stderr());
+    q(status=1)
+}
+if (is.null(opt$out))
+{
+    write("-o/--outfile missing.", stderr());
+    q(status=1)
+}
+
+# ----------------------------------------------------------------------------
+# Compute metrics.
+# ----------------------------------------------------------------------------
+
+# Read input files.x
+roi = readROI(opt$infile)
+
+# Compute some basic statistics.
+roi$min = as.numeric(lapply(roi$counts, min))
+roi$median =  as.numeric(lapply(roi$counts, median))
+roi$mean = as.numeric(lapply(roi$counts, mean))
+roi$quantile75 = as.numeric(lapply(roi$counts, quantile, probs=0.75))
+roi$quantile95 = as.numeric(lapply(roi$counts, quantile, probs=0.95))
+
+# Compute area under curve metric.
+aoc <- function(v){
+  m = max(v)
+  return (sum(v/m)/length(v))
+}
+roi$aoc = as.numeric(lapply(roi$counts, aoc))
+
+# Compute xreaXX metric.
+xreaXX <- function(v){
+  vo=order(v)
+  mx = max(v)
+  mi = min(v)
+  if((mx-mi)==0) return(0)
+  vn=(v-mi)/(mx-mi)
+  x=c(1:length(v))/length(v)
+  vno=vn[vo]
+  return(sum(x-vno)/length(v))
+}
+roi$xreaXX  =  as.numeric(lapply(roi$counts, xreaXX))
+
+# Compute r3linear metric.
+roi$r3linear = as.numeric(lapply(roi$counts, k3.linear))
+
+# Compute distMax3p metric.
+distMax3p <- function(v)
+{
+   return(length(v)-which.max(rev(v))+1)
+}
+roi$distMax3p = as.numeric(lapply(roi$count, distMax3p))
+
+# Compute distMax5p metric.
+distMax5p <- function(v)
+{
+   return (which.max(v))
+}
+roi$distMax5p = as.numeric(lapply(roi$count, distMax5p))
+
+# Write output.
+writeROI(roi, opt$outfile)
+
+# Quit with success code.
+q(status=0)
\ No newline at end of file