diff manhattan.R @ 0:c4bf5e913c2e draft default tip

"planemo upload commit b1883bac95e73fc6ffe2a36db3115ad5e5a1eba4"
author iuc
date Fri, 11 Oct 2019 17:31:09 -0400
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/manhattan.R	Fri Oct 11 17:31:09 2019 -0400
@@ -0,0 +1,60 @@
+options( show.error.messages=F, error = function () { cat( geterrmessage(), file=stderr() ); q( "no", 1, F ) } )
+
+# we need that to not crash galaxy with an UTF8 error on German LC settings.
+loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
+
+suppressPackageStartupMessages({
+    library(GWASTools)
+    library(optparse)
+})
+option_list <- list(
+    make_option(c("-f", "--file"), type="character", help="RInput GWAS file"),
+    make_option("--pval", type="integer", help="Pvalue column"),
+    make_option("--chromosome", type="integer", help="Chromosome column"),
+    make_option("--ymin", type="double", help="Min y value"),
+    make_option("--ymax", type="double", help="Max y value"),
+    make_option("--trunc", help="Show truncation lines", action="store_true"),
+    make_option("--sig", type="double", help="Significance level for lines"),
+    make_option("--thin", type="double", help="Thinning value", action="store_true", dest="thin"),
+    make_option("--ppb", type="integer", help="Points per bin, if thinning value is specified", action="store_true"))
+args <- parse_args(OptionParser(option_list=option_list))
+file <- args$file
+data <-  read.table(args$file, header=TRUE)
+pval <-  data[,args$pval]
+chromosome <-  data[,args$chromosome]
+if(!is.null(args$ymin) & !is.null(args$ymax)){
+    ylimit <-  c(args$ymin,args$ymax)
+}else if(xor(!is.null(args$ymin), !is.null(args$ymax))){
+    print("If specifying range, both ymin and ymax must be set")
+    ylimit <- NULL
+}else{
+    ylimit <- NULL
+}
+if(is.null(args$trunc)){
+    trunc <- FALSE
+}else{
+    trunc <- TRUE
+}
+if(!is.null(args$sig)){
+    sig <-  args$sig
+}else{
+    sig <- NULL
+}
+if(!is.null(args$thin)){
+    thin = args$thin
+    if(thin == 0){
+        thin = NULL 
+    }
+}else{
+    thin <- FALSE
+}
+if(!is.null(thin) & !is.null(args$ppb)){
+    ppb = args$ppb
+}
+pdf("manhattan.pdf")
+if(isFALSE(thin)){
+    manhattanPlot(pval, chromosome, ylim=ylimit, trunc.lines=trunc, signif=sig)
+}else{
+    manhattanPlot(pval, chromosome, ylim=ylimit, trunc.lines=trunc, signif=sig, thinThreshold = thin, pointsPerBin = ppb)   
+}
+invisible(dev.off())