Repository 'probmetab'
hg clone https://toolshed.g2.bx.psu.edu/repos/mmonsoor/probmetab

Changeset 0:e13ec2c3fabe (2016-07-04)
Next changeset 1:c4bd0c40eb3b (2016-07-04)
Commit message:
planemo upload commit 25fd6a739741295e3f434e0be0286dee61e06825
added:
ProbMetab.xml
README.txt
export.class.table-color-graph.R
lib.r
planemo_test.sh
probmetab.r
static/images/probmetab_workflow.png
test-data/faahKO_reduce.zip
test-data/faahOK.xset.group.retcor.group.fillPeaks.annotate.negative.Rdata
b
diff -r 000000000000 -r e13ec2c3fabe ProbMetab.xml
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/ProbMetab.xml Mon Jul 04 04:29:25 2016 -0400
[
b'@@ -0,0 +1,428 @@\n+<tool id="Probmetab" name="ProbMetab Tool" version="1.0.1">\n+\n+    <description>Wrapper function for ProbMetab R package.</description>\n+\n+    <requirements>\n+\t    <requirement type="package" version="0.4_1">r-snow</requirement>\n+        <requirement type="package" version="1.0">r-probmetab</requirement>\n+        <requirement type="package" version="1.1_4">r-batch</requirement>\n+    </requirements>\n+\n+    <stdio>\n+        <exit_code range="1:" level="fatal" />\n+    </stdio>\n+\n+    <command>\n+        LANG=C Rscript $__tool_directory__/probmetab.r \n+        #if $acquisition_options.mode == "one":\n+            mode_acquisition $acquisition_options.mode xa $acquisition_options.xa\n+            ##if $acquisition_options.xsetnofill_options.option == "show":\n+\t            ##xsetnofill $acquisition_options.xsetnofill_options.xsetnofill\n+     \t    ##end if\t\n+        #else\n+            mode_acquisition $acquisition_options.mode inputs_mode $acquisition_options.input_mode.option\n+            #if $acquisition_options.input_mode.option== "two":\n+\n+\t            image_pos $acquisition_options.input_mode.image_pos image_neg $acquisition_options.input_mode.image_neg\n+\t            ##if $acquisition_options.input_mode.xsetnofill_options.option == "show":\n+\t\t            ##xsetPnofill $acquisition_options.input_mode.xsetnofill_options.xsetPnofill xsetNnofill $acquisition_options.input_mode.xsetnofill_options.xsetNnofill\n+         \t\t##end if\t\t\n+            ##else\n+\t            ##image_combinexsannos $acquisition_options.input_mode.image_combinexsannos image_pos $acquisition_options.input_mode.image_pos\n+            #end if\n+\n+        #end if\n+\n+        #if $option_toexclude.option == "show":\n+        \ttoexclude $option_toexclude.toexclude\n+        #end if\n+        allowMiss $allowMiss html $html kegg_db $kegg_db ppm_tol $ppm_tol\n+        opt $opt corths $corths corprob $corprob pcorprob $pcorprob prob $prob\n+        \n+        #if $zip_file:\n+            zipfile $zip_file\n+        #end if\n+\n+    </command>\n+\n+    <inputs>\n+\n+        <conditional name="acquisition_options">\n+\t        <param name="mode" type="select" label="Choose your acquisition mode" >\n+\t\t        <option value="one" selected="true"  >One acquisition charge mode</option>\n+            \t<option value="two"  >Two acquisition charge mode (positif and negatif)</option>\n+\t        </param>\n+\n+\t        <!-- One acquisition mode-->\n+\t        <when value="one">\n+\t\t        <param name="xa" type="data" label="Annotate RData" format="rdata.camera.positive,rdata.camera.negative,rdata" help="Output file from annotate step " />\n+\t\t        <!--\n+\t\t        <conditional name="xsetnofill_options">\n+            \t\t\t<param name="option" type="select" label="RData group step" help="xcmsSet xcms object after missing data replacement, to retrieve SNR to isotopic peaks." >\n+             \t\t\t\t<option value="show">show</option>\n+             \t\t\t\t<option value="hide" selected="true">hide</option>\n+             \t\t\t</param>\n+             \t\t\t<when value="show">\n+\t          \t\t\t<param name="xsetnofill" type="data" label="Positive or Negative RData from group step before fillpeaks " format="rdata" help=" output from group step" />\n+             \t\t\t</when>\n+\t\t\n+          \t\t\t</conditional>\n+\t\t        -->\n+\t        </when>\n+\t        <!-- Two acquisition modes-->\n+\t        <when value="two">\n+\n+\t\t\n+\t\t        <conditional name="input_mode">\n+\t\t\t        <param name="option" type="select" label="Choose your input type method:" >\n+\t\t\t\t        <!-- Bug combinexsannos TODO <option value="one">Input from combinexsAnnos step</option> -->\n+                \t\t<option value="two" selected="true">Rdata inputs from annotate</option>\n+\t\t\t        </param>\n+\t\t\t        <!--\n+\t\t\t        <when value="one">\n+\t\t\t\t        <param name="image_combinexsannos" type="data" label="RData output from combinexsAnnos step" format="rdata" help="output file from combinexAnnos step " /> \n+\t\t\t\t        <param name="image_pos" type="data" lab'..b'--------+------------+\n+\n+\n+----------\n+Parameters\n+----------\n+\n+**Allow Miss**\n+\n+\n+Optionally retrieves peaks with no evidence of adduct or isotope and annotate them as single charged molecules [M+/-H].\n+\n+\n+**polarity**\n+\n+\n+Acquisition mode of mass spectrometer.\n+\n+\n+**Exclude samples**\n+\n+\n+Samples to be excluded of peak counting to non-annotated peak selection.\n+\n+**Calculate**\n+\n+**intervals** \n+A vector of SNR numerical intervals, to which different carbon offset should be added to predicted C-number.\n+\n+**offset**\n+\n+A vector of empirically estimated carbon offset to be added to predicted C-number.\n+\n+**massWeigth**\n+Is the contribution parameter of the probabilistic model.\n+\n+**likelihood**\n+\n+Which noise model to use, "erfc" to complementary error function, or "gaussian" to standard gaussian with two sd corresponding to the given p.p.m precision.\n+\n+**precision** \n+\n+Equipment mass accuracy, usually the same used in exact mass search.\n+\n+**KEGG database**\n+\n+\n+Select if you want to search on all KEGG organisms or multiple organisms (id1,id2,id3,...).By default,the value is KEGG which means searching on all KEGG organism. The list of KEGG IDs are available at "http://rest.kegg.jp/list/organism".\n+\n+**ppm.tol**\n+\n+\n+Parts per million mass tolerance allowed in the mass search (create.reactionMfunction).\n+\n+**HTML**\n+\n+\n+Logical, check if you want to generate a HTML ProbMetab report.This parameter uses the raw data to plot EICs and may be time consuming.\n+\n+**opt**\n+\n+\n+(reac2cor function) Correlation option, cor for correlation, and pcor for partial correlation.\n+\n+**corprob**\n+\n+(reac2cor function) Probability that the correlation is considered significant.\n+\n+\n+**pcorprob**\n+\n+\n+(reac2cor function) Probability that the partial correlation is considered significant.\n+\n+**corths**\n+\n+\n+(reac2cor function) Correlation intensity threshold.\n+\n+**prob**\n+\n+(export.class.table). How to calculate the probability to attribute a mass to a compound. Default is "count", which divide the number of times each identity was was attributed by the number of samples. Optionally the user could choose to use the mean of the probabilities of the identity, "mean".\n+\n+\n+**organismIdorganismId**\n+\n+(create.pathway.node.attributes function) KEGG organism id (http://www.kegg.jp/kegg/catalog/org_list.html) to filter possibibly pathwyas for known pathways for that organism. Only works for KEGG database for now. Default isNULL (all KEGG organisms).\n+\n+\n+\n+------------\n+Output files\n+------------\n+\n+Probmetab.RData\n+    | Rdata file, that be used outside Galaxy in R.\n+\n+Probmetab.Analysis_Report_htmlSelect if you want to search on all KEGG organisms or multiple organisms (id1,id2,id3,...).By default,the value is KEGG which means searching on all KEGG organism. The list of KEGG IDs are available at http://rest.kegg.jp/list/organism\n+    |  A list with a matrix "classTable" with attributions and probabilities and indexes of selected masses from xcms peak table (HTML format).\n+\n+Probmetab.CytoScape_output_Attribute_List.tsv\n+    |  A list with a matrix "classTable" with attributions and probabilities and indexes of selected masses from xcms peak table, that can be used as Attribute table list in CytoScape (tsv format).\n+\n+Probmetab.Analysis_Report_EICs_plots\n+    | Zip file containing the EIC plots (PNG format) of the metabolites that are listed in the HTML or tsv report.\n+\n+Probmetab.CytoScape_output.sif\n+    | Sif format file that can be used in CytoScape to visualize the network.\n+\n+\n+---------------------------------------------------\n+\n+Changelog/News\n+--------------\n+\n+**Version 1.0.1 - 16/05/2016**\n+\n+- TEST: refactoring to pass planemo test using conda dependencies\n+\n+\n+**Version 1.0.0 - 10/06/2015**\n+\n+- NEW: ProbMetab first version\n+\n+ \n+    </help>\n+    <citations>\n+        <citation type="doi">10.1093/bioinformatics/btu019</citation>\n+        <citation type="doi">10.1093/bioinformatics/btu813</citation>\n+    </citations>\n+</tool>\n+ \n+ \n'
b
diff -r 000000000000 -r e13ec2c3fabe README.txt
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/README.txt Mon Jul 04 04:29:25 2016 -0400
b
@@ -0,0 +1,19 @@
+Changelog/News
+--------------
+
+**Version 1.0.1 - 16/05/2016**
+
+- TEST: refactoring to pass planemo test using conda dependencies
+
+
+**Version 1.0.0 - 10/06/2015**
+
+- NEW: ProbMetab first version
+
+
+Test Status
+-----------
+
+Planemo test using conda: passed
+
+Planemo shed_test : passed
b
diff -r 000000000000 -r e13ec2c3fabe export.class.table-color-graph.R
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/export.class.table-color-graph.R Mon Jul 04 04:29:25 2016 -0400
[
b'@@ -0,0 +1,455 @@\n+#\' export.class.table\n+#\'\n+#\' Builds a matrix with the probability for all mass to candidate compounds\n+#\' assignments, by averaging the number of assignments obtained by the gibbs sampler algorithm\n+#\' or ordering the compound candidates with the likelihood matrix. \n+#\'\n+#\' @param gibbsL a list of attributions and probabilities from gibbs.samp function. \n+#\' @param reactionM data.frame with compound annotation information. \n+#\' @param molIon non redundant ion annotation. \n+#\' @param probM optionally to gibbsL, a matrix of likelihoods. \n+#\' @param html logical, indicating whether a html file should be generated. This parameter uses the raw data to plot EICs and may be time consuming.\n+#\' @param filename html file name, the default is "test". \n+#\' @param burnIn how many samples of the gibbs sampler should be discarded.\n+#\' @param linkPattern which pattern should be linked to compound id, for now we have \n+#\' implemented "kegg", "pubchem" and "chebi" patterns.\n+#\' @param m.test statistical test to compare mean differences. This option\n+#\' is only available to single acquisition mode analysis, with options\n+#\' "t.test" and "anova". \n+#\' @param class1 if the m.test is "t.test" first class to compare in the test,\n+#\' according with xcmsSet phenoData.\n+#\' @param class2 if the m.test is "t.test" second class to compare in the test,\n+#\' according with xcmsSet phenoData.\n+#\' @param norm logical, if TRUE performs median normalization from (Anal. Chem. 2011, 83, 5864-5872).\n+#\' @param DB data.frame table used to search compounds, with the field name to be incorporated in the html table. \n+#\' @param prob how to calculate the probability to attribute a mass to a compound. \n+#\' Default is "count", which divide the number of times each identity was \n+#\' was attributed by the number of samples. Optionally the user could\n+#\' choose to use the mean of the probabilities of the identity, "mean". \n+#\' @return A list with a matrix "classTable" with attributions and probabilities and\n+#\' indexes of selected masses from xcms peak table. \n+#\'\n+#\' @export\n+\n+export.class.table <- function(gibbsL=NULL, reactionM, molIon=NULL, probM=NULL, html=FALSE, filename="test", burnIn=3000, linkPattern="kegg", m.test="none", class1=NULL, class2=NULL, norm=FALSE, DB, prob="count", colorplot=FALSE, addLink=NULL) {       \n+\n+\tplotEIC <- function (xcmsObject, figidx, pngidx, colorplot, mode=NULL) {\n+\t\tdir.create(paste(filename,"_fig",sep=""))\n+\t\tgt<-groups(xcmsObject)\n+\t\tif(colorplot==TRUE){\n+\t\t    \tgt2 <- gt[figidx,]\n+\t\t\trgt <- gt2[,c("rtmin","rtmax")] \n+\t\t\trgt[,1] <- rgt[,1]-100 \n+\t\t\trgt[,2] <- rgt[,2]+100 \n+\t\t\t#require(doMC) \n+\t\t\t#registerDoMC() \n+\t\t\t#system.time(\n+\t\t\t#foreach(i=1:nrow(gt2)) %dopar% {\n+\t\t\tfor(i in 1:nrow(gt2)){\n+\t\t\tgroupidx1 <- which(gt[,"rtmed"] > rgt[i,1] & gt[,"rtmed"] < rgt[i,2] & gt[,"mzmed"]> gt2[i,"mzmed"]-0.1 & gt[,"mzmed"]< gt2[i,"mzmed"]+0.1)\n+\t\t\teiccor <- getEIC(xcmsObject, groupidx = groupidx1) \n+\t\t\tpng(paste(filename, "_fig/", sprintf("%003d", i), ".png", sep="")) \n+\t\t\tplot(eiccor, xcmsObject, groupidx = 1) \n+\t\t\tdev.off()\n+\t\t\t}\n+\t\t} else {\n+\t\t\tgt <- gt[figidx,]\n+\t\t\trgt <- gt[,c("rtmin","rtmax")] \n+\t\t\trgt[,1] <- rgt[,1]-100 \n+\t\t\trgt[,2] <- rgt[,2]+100 \n+\n+\t\t\teics <- getEIC(xcmsObject, mzrange=gt, rtrange =rgt, groupidx = 1:nrow(gt)) \n+\t\t\tpng(file.path(paste(filename, "_fig/%003d.png", sep="")), height=768, width=1024) \n+\t\t\t#png(file.path(paste(filename, "_fig/", pngidx, sep="")), h=768, w=1024) \n+\t\t\tplot(eics, xcmsObject)\n+\t\t\tdev.off()\n+\t\t}\n+\t\t\tif(!is.null(mode)) {\n+\t\t\t\tpngs <- dir(paste(filename, "_fig/", sep=""))\n+\t\t\t\tif(length(grep("pos|neg" , pngs))) pngs <- pngs[-grep("pos|neg" , pngs)] \n+\t\t\t\topng <- as.numeric(sub(".png","", pngs)) \n+\t\t\t\tpngs <- pngs[order(opng)]\n+\t\t\t\tname1 <- paste(filename, "_fig/",pngs, sep="") \n+\t\t\t\tname2 <- paste(filename, "_fig/",pngidx, mode, ".png", sep="")\n+\t\t\t\tfor(i in 1:length(name1)) file.rename(name1[i], name2[i])\n+\t\t\t}\n+\n+\t}\n+    allion <- molIon$molIon[molIon$molIon[,"isotope"]==0'..b's & Samples]\n+\t   medSam <- apply(mat, 1, median)\n+\t   medSam[which(medSam==0)] <- 0.0001\n+\t   mat <- apply(mat, 2, function(mat, medSam){\n+\t      medFDiSmpl <- mat/medSam\n+\t      vec<-mat/median(medFDiSmpl)\n+\t      return(vec)\n+\t   }, medSam)\n+\t   return (mat)\n+\t}\n+\t# this piece of code was copied from xcms\n+\tpval <- function(X, classlabel, teststat) {\n+\n+\t    n1 <- rowSums(!is.na(X[,classlabel == 0]))\n+\t    n2 <- rowSums(!is.na(X[,classlabel == 1]))\n+\t    A <- apply(X[,classlabel == 0], 1, sd, na.rm=TRUE)^2/n1 ## sd(t(X[,classlabel == 0]), na.rm = TRUE)^2/n1\n+\t    B <- apply(X[,classlabel == 1], 1, sd, na.rm=TRUE)^2/n2 ## sd(t(X[,classlabel == 1]), na.rm = TRUE)^2/n2\n+\t    df <- (A+B)^2/(A^2/(n1-1)+B^2/(n2-1))\n+\n+\t    pvalue <- 2 * (1 - pt(abs(teststat), df))\n+\t    invisible(pvalue)\n+\t}\n+\n+\tc1 <- grep(class1, molIon$cameraobj@xcmsSet@phenoData[,1]) \n+\tc2 <- grep(class2, molIon$cameraobj@xcmsSet@phenoData[,1]) \n+\ttestclab <- c(rep(0, length(c1)), rep(1, length(c2)))\n+\ttestval <- groupval(molIon$cameraobj@xcmsSet, "medret", "into") \n+\tif(norm) testval <- normalize.medFC(testval) \t\n+\ttstat <- mt.teststat(testval, testclab)\n+        pvalue <- pval(testval, testclab, tstat)\n+\n+#\n+#\t\trport <- diffreport(molIon$cameraobj@xcmsSet, class1=class1, class2= class2, sortpval=FALSE) \n+#\t\tans[ans[,6]!="",5] <- rport[figidx, "pvalue"]\n+\t\tans[ans[,6]!="",5] <- pvalue[figidx]\n+\t\ttestname <- "t.test p-value"\n+\t}\n+\tif(m.test=="anova"){\n+\t\tclass <- molIon$cameraobj@xcmsSet@phenoData \n+\t\tgetPvalue <- function(dataidx) {\n+\t\t\taov.data <- data.frame(resp=as.numeric(peaklist[dataidx,coords]), class=class) \n+\t\t\tanova(aov(resp~class, aov.data))$Pr[1] \n+\t\t}\n+\t\ttestname <- "anova p-value"\n+\t\tans[ans[,6]!="",5] <- sapply(figidx, getPvalue)\n+\t}\n+\n+\theader <- matrix(c("Proposed Mass","Most probable Compound","Probability","Entropy", testname,"EIC-plot", "Ion annotation",coords), nrow=1 , ncol=7+length(coords) ) \n+\tans <- rbind(header, ans)\n+\n+\n+\t# additional field\n+\t# ans <- cbind(ans[,1:2], ans[,2], ans[,3:ncol(ans)]) \n+\t#ans[ans[,3]!="unknown",][-1,3] <- as.character(DB$sbml.id[nid]) \n+\n+\tif(html) {\n+\t\t#require(hwriter)\n+\t\tansb <- ans\n+\t\tans[ans[,2]!="unknown",][-1,2] <- as.character(DB$name[nid]) \n+\t\tif(linkPattern=="pubchem") ansb <- ans\n+\t\t\n+\t\thyper=matrix(paste(linkURL, ansb[-1,2], sep=""),ncol=1 )\n+\t\tif(!is.null(molIon$cameraobj)) {\t\n+\t\t\thyper1=matrix(paste(fig, ans[-1,6],".png", sep=""),ncol=1 ) \t\n+\t\t}\n+\t\telse {\n+\t\t\thyper1 <- ans[-1,6]\n+\t\t\thyper1[ans[-1,6]!=""][allidx[,2]!="neg"] <- paste(hyper1[ans[-1,6]!=""][allidx[,2]!="neg"], "pos", sep="")\n+\t\t\thyper1[ans[-1,6]!=""][allidx[,2]=="neg"] <- paste(hyper1[ans[-1,6]!=""][allidx[,2]=="neg"], "neg", sep="")\n+\t\t\thyper1=matrix(paste(fig, hyper1,".png", sep=""),ncol=1 ) \t\n+\t\t}\t\n+\t\tp=openPage(paste(filename,".html",sep=""))\n+\t\tans2 <- ans[,1:7]\n+\t\tlink <- cbind(matrix(NA,nrow(ans2),1),rbind(NA,hyper),matrix(NA,nrow(ans2),3),rbind(NA,hyper1),matrix(NA,nrow(ans2),1))\n+\t\t# This block is responsible to add as many columns as the user\n+\t\t# wants\n+\t\tif(!is.null(addLink)){ \n+\t\t\tfor(l in 1:length(addLink$link)) {\n+\t\t\t\tlink <- cbind(link, rbind(NA, addLink[[l]]))\n+\t\t\t}\n+\t\t\tfor(a in 1:length(addLink$ans)) {\n+\t\t\t\tans2 <- cbind(ans2,addLink$ans[[a]]) \n+\t\t\t\t#colnames(ans2)[7+a] <- addLink$ans[[a]][1]\n+\t\t\t}\n+\t\t}\n+\t\thwrite(ans2, p,row.bgcolor=\'#ffdc98\', link=link  )\n+\t\tclosePage(p)\n+\t\tif(!is.null(molIon$cameraobj)) {\t\n+\t\t\tplotEIC(molIon$cameraobj@xcmsSet, figidx, ans[ans[,6]!="",6][-1], colorplot=colorplot) \n+\t\t}\t\n+\t\telse {\n+\t\t\t dataidxP <- as.numeric(allidx[allidx[,2]!="neg",1]) \n+\t\t\t pngidxP <- allidx[allidx[,2]!="neg",3]\n+\t\t\tplotEIC(molIon$pos@xcmsSet, dataidxP, pngidxP, "pos", colorplot=colorplot) \n+\t\t\t dataidxN <- as.numeric(allidx[allidx[,2]=="neg",1]) \n+\t\t\t pngidxN <- allidx[allidx[,2]=="neg",3]\n+\t\t\tplotEIC(molIon$neg@xcmsSet, dataidxN, pngidxN, "neg", colorplot=colorplot) \n+\t\t}\n+\t\t\n+\t}\n+\telse {\n+\t\tansb <- ans\n+\t}\n+\tcolnames(ansb) <- ansb[1,]\n+\tansb <- ansb[-1,]\n+\treturn(list(classTable=ansb,  figidx=figidx))\n+} \n'
b
diff -r 000000000000 -r e13ec2c3fabe lib.r
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/lib.r Mon Jul 04 04:29:25 2016 -0400
[
b'@@ -0,0 +1,323 @@\n+# lib.r ProbMetab version="1.0.0"\r\n+# Author: Misharl Monsoor ABIMS TEAM mmonsoor@sb-roscoff.fr\r\n+# Contributors: Yann Guitton and Jean-francois Martin\r\n+\r\n+\r\n+##Main probmetab function launch by the Galaxy ProbMetab wrapper\r\n+probmetab = function(xa, xaP, xaN, variableMetadata, variableMetadataP, variableMetadataN, listArguments){\r\n+\t##ONE MODE ACQUISITION##\r\n+\tif(listArguments[["mode_acquisition"]]=="one") {\r\n+\t\tcomb=NULL\n+\r\n+\t\t#Get the polarity from xa object\r\n+\t\tpolarity=xa@polarity\r\n+\t\t#SNR option\r\n+\t\tif ("xsetnofill" %in% names(listArguments)) {\r\n+\t\t\tload(listArguments[["xsetnofill"]])\r\n+\t\t\txsetnofill=xset\r\n+\t\t}\r\n+\t\telse{\r\n+\t\t\txsetnofill=NULL\r\n+\t\t}\r\n+\t\t#Exclude samples\r\n+\t\tif ("toexclude" %in% names(listArguments)) {\r\n+\t\t\ttoexclude=listArguments[["toexclude"]]\r\n+\t\t}\r\n+\t\telse {\r\n+\t\t\ttoexclude=NULL\r\n+\t\t}\r\n+\t\tionAnnot=get.annot(xa, polarity=polarity, allowMiss=listArguments[["allowMiss"]],xset=xsetnofill,toexclude=toexclude)\r\n+\t\tcomb=NULL\r\n+\t}\r\n+\r\n+\t##TWO MODES ACQUISITION##\r\n+\t#Mode annotatediffreport\r\n+\telse if(listArguments[["inputs_mode"]]=="two"){\r\n+\t\t##Prepare the objects that will be used for the get.annot function\r\n+\t\tcomb=1\r\n+\t\t\r\n+\r\n+\t\txsetPnofill=NULL\r\n+\t\txsetNnofill=NULL\r\n+\t\t# TODO: a reactiver\t\t\n+\t\t#if ("xsetPnofill" %in% names(listArguments)) {\r\n+\t\t#\tload(listArguments[["xsetPnofill"]])\r\n+\t\t#\txsetPnofill=xset\r\n+\t\t#}\r\n+\t\t#if ("xsetNnofill" %in% names(listArguments)) {\r\n+\t\t#\tload(listArguments[["xsetNnofill"]])\r\n+\t\t#\txsetNnofill=xset\r\n+\t\t#}\r\n+\t\t# include CAMERA non-annotated compounds, and snr retrieval \r\n+\t\t# comb 2+ - used on Table 1\t\r\n+\t\tionAnnotP2plus = get.annot(axP, allowMiss=listArguments[["allowMiss"]], xset=xsetPnofill,toexclude=listArguments[["toexclude"]]) \r\n+\t\tionAnnotN2plus = get.annot(axN, polarity="negative", allowMiss=listArguments[["allowMiss"]], xset=xsetNnofill,toexclude=listArguments[["toexclude"]])\r\n+\t\tionAnnot = combineMolIon(ionAnnotP2plus, ionAnnotN2plus)\r\n+\t\tprint(sum(ionAnnot$molIon[,3]==1))\r\n+\t\tprint(sum(ionAnnot$molIon[,3]==0))\r\n+\t\twrite.table(ionAnnot[1], sep="\\t", quote=FALSE, row.names=FALSE, file="CombineMolIon.tsv")\r\n+\t\t#Merge variableMetadata Negative and positive acquisitions mode\r\n+\t\t\r\n+\r\n+\t\t#Mode combinexsannos TODO bug avec tableau issus de combinexsannos\r\n+\t\t#else {\r\n+\t\t\t#load(listArguments[["image_combinexsannos"]])\r\n+\t\t\t#image_combinexsannos=cAnnot\r\n+\t\t\t##Prepare the objects that will be used for the combineMolIon function\r\n+\t\t\t#load(listArguments[["image_pos"]])\r\n+\t\t\t#image_pos=xa\r\n+\t\t\t#ionAnnot=combineMolIon(peaklist=cAnnot, cameraobj=image_pos, polarity="pos")\r\n+\t\t#}\r\n+\t\t\r\n+\t}\r\n+\r\n+\t##DATABASE MATCHING##\r\n+\tif (listArguments[["kegg_db"]]=="KEGG"){\t\t\r\n+\t\tDB=build.database.kegg(orgID = NULL)\r\n+\t}\r\n+\telse{\t\r\n+\t\ttable_list <<- NULL\r\n+\t\tids=strsplit(listArguments[["kegg_db"]],",")\r\n+\t\tids=ids[[1]]\r\n+\t\tif(length(ids)>1){\r\n+\t\t\tfor(i in 1:length(ids)){\r\n+\t\t\t\t table_list[[i]] <- build.database.kegg(ids[i])\r\n+\t\t\t}\r\n+\t\t\tdb_table=do.call("rbind",table_list)\r\n+\t\t\tDB=unique(db_table)\r\n+\t\t}\r\n+\t\telse{\r\n+\t\t\tDB=build.database.kegg(listArguments[["kegg_db"]])\r\n+\t\t}\r\n+\t}\t\r\n+\t#Matching des mass exactes mesurees avec les masses des compounds KEGG (pas M+H ou M-H)\r\n+\treactionM = create.reactionM(DB, ionAnnot, ppm.tol=listArguments[["ppm_tol"]])\r\n+\t##PROBABILITY RANKING##\r\n+\t# number of masses with candidates inside the fixed mass window\r\n+\t# and masses with more than one candidate\r\n+\tlength(unique(reactionM[reactionM[,"id"]!="unknown",1])) \r\n+\tsum(table(reactionM[reactionM[,"id"]!="unknown",1])>1)\r\n+\t#if (listArguments[["useIso"]]){\r\n+\t\t#BUG TODO\r\n+\t\t# Calculate the ratio between observed and theoretical isotopic patterns.\r\n+\t\t# If you don\'t have an assessment of carbon offset to carbon number prediction \r\n+\t\t# skip this step and use the reactionM as input to weigthM function. \r\n+\t\t#isoPatt < incorporate.isotopes(comb2plus, reactionM, , samp=12:23, DB=DB)  \r\n+\t\t#  calculate   the   likelihood   of   each   mass   to   compound   assignment   using   mass   accuracy,and iso'..b'+",anstr,fixed=T)==T) || (grepl("]3+",anstr,fixed=T)==T)) { esi <- "p"}\r\n+      else \r\n+        if ((grepl("]-",anstr,fixed=T)==T) || (grepl("]2-",anstr,fixed=T)==T) || (grepl("]3-",anstr,fixed=T)==T)) { esi <- "n"}\r\n+      # cat(" ioni ",esi,"\\n")\r\n+    } else\r\n+    { esi <- "u"} \r\n+  \r\n+  return(esi)\r\n+}\r\n+\r\n+\r\n+merge_probmetab <- function(metaVar,ansConn) {\r\n+  ## Parse ProbMetab information result file and merge in variable_metaData initial file\r\n+  ##  inputs : \r\n+  ##      metaVar : data.frame of metadataVariable input of probmetab function\r\n+  ##     ansConn  : data.frame of ProbMetab result\r\n+  ## output : dataframe with Probmetab results merge with variableMetadata\r\n+  ## Constante\r\n+  ## iannot : indice de la colonne annotation dans le resultat de probMetab\r\n+  iannot <- 4\r\n+  \r\n+  ## definition of an unique identification of ions mz with 3 decimals and rt(sec) with 1 decimal to avoid \r\n+  ## duplicate ions name in the diffreport result file\r\n+  ions <- paste ("M",round(metaVar$mz,3),"T",round(metaVar$rt,1),sep="")\r\n+  metaVar <- data.frame(ions,metaVar)\r\n+  \r\n+  ###### Result data.frame from ProbMetab result list\r\n+  an_ini <- ansConn$classTable\r\n+  \r\n+  ## Suppression of rows without  mz and rt or unknown and columns of intensities\r\n+  ## COLUMNS SUBSCRIPTS HAVE TO BE CHECKED WITh DIFFERENT RESULTS FILES \r\n+  an <- an_ini[(an_ini[,2]!="unknown"),c(1,2,3,7)]\r\n+  ## initialisation of vectors receiving the result of the parse of the column annotation (subscrip iannot) \r\n+  mz <- rep(0,dim(an)[1])\r\n+  rt <- rep(0,dim(an)[1])\r\n+  propmz <- rep(0,dim(an)[1])\r\n+  ioni <- rep("u",dim(an)[1])\r\n+  \r\n+  ## parse the column annotation and define ionisation mode\r\n+  for (i in 1:dim(an)[1]) {\r\n+    if (an[i,1] != "") {\r\n+      info_mzrt <- unlist(strsplit(an[i,iannot],"#"))\r\n+      propmz[i] <- as.numeric(an[i,1])\r\n+      mz[i] <- as.numeric(info_mzrt[1])\r\n+      rt[i] <- as.numeric(info_mzrt[2])\r\n+      ioni[i] <- deter_ioni(info_mzrt,as.numeric(an[i,1]))\r\n+    }\r\n+    else {\r\n+      propmz[i] <- as.numeric(propmz[i-1])\r\n+      mz[i] <- as.numeric(mz[i-1])\r\n+      rt[i] <- as.numeric(rt[i-1])\r\n+      ioni[i] <- ioni[i-1]\r\n+    }\r\n+  }\r\n+  \r\n+  ## definition of an unique identification of ions : mz with 3 decimals and rt(sec) with 1 decimal\r\n+  ## The same as for the metadataVariable data.frame to match with.\r\n+  ions <- paste ("M",round(mz,3),"T",round(rt,1),sep="")\r\n+  an <- data.frame(ions,ioni,propmz,mz,rt,an)\r\n+  \r\n+  ## transposition of the different probmetab annotations which are in different rows in the initial result data.frame\r\n+  ## on only 1 row separated with a ";"\r\n+  li <- as.matrix(table(an$propmz))\r\n+  li <- data.frame(dimnames(li)[1],li)\r\n+  dimnames(li)[[2]][1] <- "propmz"\r\n+  ions   <- rep("u",dim(li)[1])\r\n+  propmz <- rep(0,dim(li)[1])\r\n+  mpc    <- rep("c",dim(li)[1])\r\n+  proba  <- rep("p",dim(li)[1])\r\n+  c <- 0\r\n+  while (c < dim(li)[1]) {\r\n+    c <- c + 1\r\n+    suban     <- an[an$propmz==li[c,1],]\r\n+    ions[c]   <- as.character(suban[1,1])\r\n+    propmz[c] <- as.numeric(suban[1,3])\r\n+    mpc[c]    <- paste(suban[,7],collapse=";")\r\n+    proba[c]  <- paste(as.character(suban[,8]),collapse=";") \r\n+  }\r\n+  \r\n+  ## Creation of the data.frame with 1 row per ions\r\n+  anc <- data.frame(ions,propmz,mpc,proba)\r\n+  anc <- anc[order(anc[,1]),]\r\n+  \r\n+  metaVarFinal <- merge(metaVar, anc, by.x=1, by.y=1, all.x=T, all.y=T)\r\n+  metaVarFinal <- metaVarFinal[,-1]\r\n+  #write.table(metaVarFinal,file="res.txt", sep="\\t", row.names=F, quote=F)\r\n+  \r\n+  return (metaVarFinal)\r\n+}\r\n+\n+# RETROCOMPATIBILITE avec ancienne version de annotate\n+getVariableMetadata = function(xa) {\n+\t# --- variableMetadata ---\n+\tpeakList=getPeaklist(xa)\n+\tpeakList=cbind(groupnames(xa@xcmsSet),peakList); colnames(peakList)[1] = c("name");\n+\tvariableMetadata=peakList[,!(colnames(peakList) %in% c(sampnames(xa@xcmsSet)))]\n+\tvariableMetadata$name= paste("M",round(variableMetadata$mz),"T",round(variableMetadata$rt),sep="")\n+\treturn (variableMetadata)\n+}\r\n'
b
diff -r 000000000000 -r e13ec2c3fabe planemo_test.sh
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/planemo_test.sh Mon Jul 04 04:29:25 2016 -0400
b
@@ -0,0 +1,5 @@
+
+# -- Use of conda dependencies
+planemo conda_init --conda_prefix /tmp/mc
+planemo conda_install --conda_prefix /tmp/mc .
+planemo test --install_galaxy --conda_prefix /tmp/mc --conda_dependency_resolution 
b
diff -r 000000000000 -r e13ec2c3fabe probmetab.r
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/probmetab.r Mon Jul 04 04:29:25 2016 -0400
[
@@ -0,0 +1,119 @@
+#!/usr/local/public/bin/Rscript --vanilla --slave --no-site-file
+# probmetab.r version="1.0.0"
+# Author: Misharl Monsoor ABIMS TEAM mmonsoor@sb-roscoff.fr
+
+
+# ----- LOG -----
+log_file=file("probmetab.log", open = "wt")
+sink(log_file)
+sink(log_file, type = "out")
+
+# ----- PACKAGE -----
+cat("\tPACKAGE INFO\n")
+pkgs=c("parallel","BiocGenerics", "Biobase", "Rcpp", "mzR", "igraph", "xcms","snow","CAMERA","batch","ProbMetab")
+for(p in pkgs) {
+ suppressWarnings( suppressPackageStartupMessages( stopifnot( library(p, quietly=TRUE, logical.return=TRUE, character.only=TRUE))))
+ cat(p,"\t",as.character(packageVersion(p)),"\n",sep="")
+}
+
+source_local <- function(fname){
+    argv <- commandArgs(trailingOnly = FALSE)
+    base_dir <- dirname(substring(argv[grep("--file=", argv)], 8))
+    source(paste(base_dir, fname, sep="/"))
+}
+cat("\n\n")
+# ----- ARGUMENTS -----
+cat("\tARGUMENTS INFO\n") 
+listArguments = parseCommandArgs(evaluate=FALSE) #interpretation of arguments given in command line as an R list of objects
+write.table(as.matrix(listArguments), col.names=F, quote=F, sep='\t')
+
+
+
+# ----- PROCESSING INFILE -----
+cat("\tINFILE PROCESSING INFO\n")
+
+# ----- INFILE PROCESSING -----
+
+if(listArguments[["mode_acquisition"]]=="one") {
+ load(listArguments[["xa"]])
+    
+    if (!is.null(listArguments[["zipfile"]])){
+        zipfile= listArguments[["zipfile"]]; listArguments[["zipfile"]]=NULL
+    }
+
+    #Unzip the chromatograms file for plotting EIC pour the HTML file
+ if(exists("zipfile"))
+ {
+ if (zipfile!=""){
+ directory=unzip(zipfile)
+ }
+ }
+ if (!exists("xa")) {
+ xa=xsAnnotate_object
+ }
+ source_local("lib.r")
+ if (!exists("variableMetadata")) variableMetadata= getVariableMetadata(xa);
+
+} else if(listArguments[["inputs_mode"]]=="two"){
+ load(listArguments[["image_pos"]])
+    
+    if (!is.null(listArguments[["zipfile"]])){
+        zipfile= listArguments[["zipfile"]]; listArguments[["zipfile"]]=NULL
+    }
+    
+ #Unzip the chromatograms file for plotting EIC pour the HTML file
+ if(exists("zipfile")) {
+ if (zipfile!=""){
+ directory=unzip(zipfile)
+ }
+ }
+ if (!exists("xa")) {
+ xa=xsAnnotate_object
+ }
+ xaP=xa
+ source_local("lib.r")
+ if (!exists("variableMetadata")) variableMetadataP= getVariableMetadata(xa)
+ else variableMetadataP=variableMetadata
+
+
+ load(listArguments[["image_neg"]])
+    
+    if (!is.null(listArguments[["zipfile"]])){
+        zipfile= listArguments[["zipfile"]]; listArguments[["zipfile"]]=NULL
+    }
+    
+ #Unzip the chromatograms file for plotting EIC pour the HTML file
+ if(exists("zipfile")) {
+
+   if (zipfile!=""){
+ directory=unzip(zipfile)
+ }
+ }
+ if (!exists("xa")) {
+ xa=xsAnnotate_object
+ }
+ xaN=xa
+ source_local("lib.r")
+
+ if (!exists("variableMetadata")) variableMetadataN= getVariableMetadata(xa)
+ else variableMetadataN=variableMetadata
+}
+
+#Import the different functions
+source_local("lib.r")
+source_local("export.class.table-color-graph.R")
+
+# ----- PROCESSING INFO -----
+cat("\tMAIN PROCESSING INFO\n")
+
+if(listArguments[["mode_acquisition"]]=="one") {
+ results=probmetab(xa=xa, variableMetadata=variableMetadata,listArguments=listArguments)
+} else if(listArguments[["inputs_mode"]]=="two"){
+ results=probmetab(xaP=xaP, xaN=xaN,variableMetadataP=variableMetadataP, variableMetadataN=variableMetadataN, listArguments=listArguments)
+}
+#delete the parameters to avoid the passage to the next tool in .RData image
+#rm(listArguments)
+cat("\tDONE\n")
+#saving R data in .Rdata file to save the variables used in the present tool
+#save.image(paste("probmetab","RData",sep="."))
+
b
diff -r 000000000000 -r e13ec2c3fabe static/images/probmetab_workflow.png
b
Binary file static/images/probmetab_workflow.png has changed
b
diff -r 000000000000 -r e13ec2c3fabe test-data/faahKO_reduce.zip
b
Binary file test-data/faahKO_reduce.zip has changed
b
diff -r 000000000000 -r e13ec2c3fabe test-data/faahOK.xset.group.retcor.group.fillPeaks.annotate.negative.Rdata
b
Binary file test-data/faahOK.xset.group.retcor.group.fillPeaks.annotate.negative.Rdata has changed