diff roi_pca.Rscript @ 0:61d9bdb6d519 draft

Uploaded

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/roi_pca.Rscript	Thu Apr 18 08:03:38 2013 -0400
@@ -0,0 +1,66 @@
+#!/usr/bin/env Rscript
+
+# Compute PCA metric for ROI files.
+#
+# USAGE: roi_pca.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 PCA.
+data=roi[,8:length(roi)-1]
+pca=princomp(data, na.action=na.exclude)
+prj=predict(pca, data)
+roi$pca1=prj[,1]
+roi$pca2=prj[,2]
+
+# Write output.
+writeROI(roi, opt$outfile)
+
+# Quit with success code.
+q(status=0)
\ No newline at end of file