Repository 'limma_voom'
hg clone https://toolshed.g2.bx.psu.edu/repos/iuc/limma_voom

Changeset 4:a61a6e62e91f (2018-04-29)
Previous changeset 3:38aab66ae5cb (2018-01-31) Next changeset 5:d8a55b5f0de0 (2018-05-05)
Commit message:
planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/limma_voom commit 6a458881c0819b75e55e64b3f494679d43bb9ee8
modified:
limma_voom.R
limma_voom.xml
test-data/anno.txt
added:
test-data/out_rscript.txt
removed:
test-data/limma-trend_Mut-WT.tsv
test-data/limma-voom_Mut-WT.tsv
test-data/limma-voom_Mut-WT_2fact.tsv
test-data/limma-voom_Mut-WT_2fact_anno.tsv
test-data/limma-voom_Mut-WT_anno.tsv
test-data/limma-voom_WT-Mut.tsv
test-data/limma-voom_WT-Mut_2fact_anno.tsv
test-data/limma-voom_normcounts.tsv
test-data/limma-voom_normcounts_anno.tsv
b
diff -r 38aab66ae5cb -r a61a6e62e91f limma_voom.R
--- a/limma_voom.R Wed Jan 31 12:45:42 2018 -0500
+++ b/limma_voom.R Sun Apr 29 17:36:42 2018 -0400
[
@@ -349,7 +349,9 @@
 data <- list()
 data$counts <- counts
 if (haveAnno) {
-  data$genes <- geneanno
+  # order annotation by genes in counts (assumes gene ids are in 1st column of geneanno)
+  annoord <- geneanno[match(row.names(counts), geneanno[,1]), ]
+  data$genes <- annoord
 } else {
   data$genes <- data.frame(GeneID=row.names(counts))
 }
@@ -431,6 +433,7 @@
     # limma-trend approach
     logCPM <- cpm(data, log=TRUE, prior.count=opt$trend)
     fit <- lmFit(logCPM, design)
+    fit$genes <- data$genes
     fit <- contrasts.fit(fit, contrasts)
     if (wantRobust) {
         fit <- eBayes(fit, trend=TRUE, robust=TRUE)
@@ -459,7 +462,7 @@
 
     # Save normalised counts (log2cpm)
     if (wantNorm) {
-        write.table(logCPM, file=normOut, row.names=TRUE, sep="\t")
+        write.table(logCPM, file=normOut, row.names=TRUE, sep="\t", quote=FALSE)
         linkData <- rbind(linkData, c((paste0(deMethod, "_", "normcounts.tsv")), (paste0(deMethod, "_", "normcounts.tsv"))))
     }
 } else {
@@ -510,7 +513,7 @@
      # Save normalised counts (log2cpm)
     if (wantNorm) {
         norm_counts <- data.frame(vData$genes, vData$E)
-        write.table(norm_counts, file=normOut, row.names=FALSE, sep="\t")
+        write.table(norm_counts, file=normOut, row.names=FALSE, sep="\t", quote=FALSE)
         linkData <- rbind(linkData, c((paste0(deMethod, "_", "normcounts.tsv")), (paste0(deMethod, "_", "normcounts.tsv"))))
     }
 
@@ -554,11 +557,7 @@
 
     # Write top expressions table
     top <- topTable(fit, coef=i, number=Inf, sort.by="P")
-    if (wantTrend) {
-        write.table(top, file=topOut[i], row.names=TRUE, sep="\t")
-    } else {
-        write.table(top, file=topOut[i], row.names=FALSE, sep="\t")
-    }
+    write.table(top, file=topOut[i], row.names=FALSE, sep="\t", quote=FALSE)
 
     linkName <- paste0(deMethod, "_", contrastData[i], ".tsv")
     linkAddr <- paste0(deMethod, "_", contrastData[i], ".tsv")
@@ -627,7 +626,7 @@
 
 cata("<body>\n")
 cata("<h3>Limma Analysis Output:</h3>\n")
-cata("PDF copies of JPEGS available in 'Plots' section.<br />\n")
+cata("Links to PDF copies of plots are in 'Plots' section below />\n")
 if (wantWeight) {
     HtmlImage(imageData$Link[1], imageData$Label[1], width=1000)
 } else {
b
diff -r 38aab66ae5cb -r a61a6e62e91f limma_voom.xml
--- a/limma_voom.xml Wed Jan 31 12:45:42 2018 -0500
+++ b/limma_voom.xml Sun Apr 29 17:36:42 2018 -0400
[
b'@@ -1,10 +1,10 @@\n-<tool id="limma_voom" name="limma" version="3.34.6.0">\n+<tool id="limma_voom" name="limma" version="3.34.9.0">\n     <description>\n         Perform differential expression with limma-voom or limma-trend\n     </description>\n \n     <requirements>\n-        <requirement type="package" version="3.34.6">bioconductor-limma</requirement>\n+        <requirement type="package" version="3.34.9">bioconductor-limma</requirement>\n         <requirement type="package" version="3.20.7">bioconductor-edger</requirement>\n         <requirement type="package" version="1.4.30">r-statmod</requirement>\n         <requirement type="package" version="0.5.0">r-scales</requirement>\n@@ -102,6 +102,10 @@\n \n &&\n cp \'$outReport.files_path\'/*.tsv output_dir/\n+\n+#if $out.rscript:\n+    && cp \'$__tool_directory__/limma_voom.R\' \'$rscript\'\n+#end if\n     ]]></command>\n \n     <inputs>\n@@ -138,12 +142,12 @@\n                     </param>\n                     <repeat name="rep_group" title="Group" min="2" default="2">\n                         <param name="groupName" type="text" label="Name"\n-                        help="Name of group that the counts files(s) belong to (e.g. WT or Mut). NOTE: Please only use letters, numbers or underscores (case sensitive).">\n+                        help="Name of group that the counts files belong to (e.g. WT or Mut). NOTE: Please only use letters, numbers or underscores (case sensitive).">\n                         <sanitizer>\n                             <valid initial="string.letters,string.digits"><add value="_" /></valid>\n                         </sanitizer>\n                         </param>\n-                        <param name="countsFile" type="data" format="tabular" multiple="true" label="Counts file(s)"/>\n+                        <param name="countsFile" type="data" format="tabular" multiple="true" label="Counts files"/>\n                     </repeat>\n                 </repeat>\n             </when>\n@@ -245,6 +249,7 @@\n                 label="Output Normalised Counts Table?"\n                 help="Output a file containing the normalised counts, these are in log2 counts per million (logCPM). Default: No">\n             </param>\n+            <param name="rscript" type="boolean" truevalue="True" falsevalue="False" checked="False" label="Output Rscript?" help="If this option is set to Yes, the Rscript used will be provided as a text file in the output. Default: No"/>\n             <param name="rdaOption" type="boolean" truevalue="1" falsevalue="0" checked="false"\n                 label="Output RData file?"\n                 help="Output all the data used by R to construct the plots and tables, can be loaded into R. A link to the RData file will be provided in the HTML report. Default: No">\n@@ -281,6 +286,9 @@\n         <collection name="outTables" type="list" label="${tool.name} on ${on_string}: Tables">\n             <discover_datasets pattern="(?P&lt;name&gt;.+)\\.tsv$" format="tabular" directory="output_dir" visible="false" />\n         </collection>\n+        <data name="rscript" format="txt" label="${tool.name} on ${on_string}: Rscript">\n+            <filter>out[\'rscript\']</filter>\n+        </data>\n     </outputs>\n \n     <tests>\n@@ -300,8 +308,18 @@\n             </repeat>\n             <param name="normalisationOption" value="TMM" />\n             <output_collection name="outTables" count="2">\n-                <element name="limma-voom_Mut-WT" ftype="tabular" file="limma-voom_Mut-WT.tsv" />\n-                <element name="limma-voom_WT-Mut" ftype="tabular" file="limma-voom_WT-Mut.tsv" />\n+                <element name="limma-voom_Mut-WT" ftype="tabular" >\n+                    <assert_contents>\n+                        <has_text_matching expression="GeneID.*logFC.*AveExpr.*t.*P.Value.*adj.P.Val.*B" />\n+                        <has_text_matching expression="11304.*0.4573" />\n+                    </assert_contents>\n+                </element>\n+                <element name="limma-voom_WT-Mut" ftype="tabular" >\n+           '..b'`_.\n \n **Counts Data:**\n-\n The counts data can either be input as separate counts files (one sample per file) or a single count matrix (one sample per column). The rows correspond to genes, and columns correspond to the counts for the samples. Values must be tab separated, with the first row containing the sample/column labels and the first column containing the row/gene labels. Gene identifiers can be of any type but must be unique and not repeated within a counts file.\n \n Example - **Separate Count Files**:\n@@ -520,19 +614,19 @@\n **Gene Annotations:**\n Optional input for gene annotations, this can contain more\n information about the genes than just an ID number. The annotations will\n-be available in the differential expression results table and the optional normalised counts table.\n+be available in the differential expression results table and the optional normalised counts table. The file must contain a header row and have the gene IDs in the first column. The number of rows should match that of the counts files, add NA for any gene IDs with no annotation. The Galaxy tool **annotateMyIDs** can be used to obtain annotations for human, mouse, fly and zebrafish.\n \n Example:\n \n     ==========  ==========  ===================================================\n     **GeneID**  **Symbol**  **GeneName**\n     ----------  ----------  ---------------------------------------------------\n-    1287        Pzp         pregnancy zone protein\n-    1298        Aanat       arylalkylamine N-acetyltransferase\n-    1302        Aatk        apoptosis-associated tyrosine kinase\n-    1303        Abca1       ATP-binding cassette, sub-family A (ABC1), member 1\n-    1304        Abca4       ATP-binding cassette, sub-family A (ABC1), member 4\n-    1305        Abca2       ATP-binding cassette, sub-family A (ABC1), member 2\n+    11287       Pzp         pregnancy zone protein\n+    11298       Aanat       arylalkylamine N-acetyltransferase\n+    11302       Aatk        apoptosis-associated tyrosine kinase\n+    11303       Abca1       ATP-binding cassette, sub-family A (ABC1), member 1\n+    11304       Abca4       ATP-binding cassette, sub-family A (ABC1), member 4\n+    11305       Abca2       ATP-binding cassette, sub-family A (ABC1), member 2\n     ==========  ==========  ===================================================\n \n **Factor Information:**\n@@ -556,24 +650,6 @@\n *Groups:* The names of the groups for the factor. These must be entered in the same order as the samples (to which the groups correspond) are listed in the columns of the counts matrix. Spaces must not be used and if entered into the tool form above, the values should be separated by commas.\n \n \n-**Gene Annotations:**\n-Optional input for gene annotations, this can contain more\n-information about the genes than just an ID number. The annotations will\n-be available in the differential expression results table and the optional normalised counts table.\n-\n-Example:\n-\n-    ==========  ==========  ===================================================\n-    **GeneID**  **Symbol**  **GeneName**\n-    ----------  ----------  ---------------------------------------------------\n-    1287        Pzp         pregnancy zone protein\n-    1298        Aanat       arylalkylamine N-acetyltransferase\n-    1302        Aatk        apoptosis-associated tyrosine kinase\n-    1303        Abca1       ATP-binding cassette, sub-family A (ABC1), member 1\n-    1304        Abca4       ATP-binding cassette, sub-family A (ABC1), member 4\n-    1305        Abca2       ATP-binding cassette, sub-family A (ABC1), member 2\n-    ==========  ==========  ===================================================\n-\n **Contrasts of Interest:**\n The contrasts you wish to make between levels.\n A common contrast would be a simple difference between two levels: "Mut-WT"\n@@ -675,6 +751,7 @@\n Optionally, under **Output Options** you can choose to output\n \n     * a normalised counts table\n+    * the R script used by this tool\n     * an RData file\n \n -----\n'
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/anno.txt
--- a/test-data/anno.txt Wed Jan 31 12:45:42 2018 -0500
+++ b/test-data/anno.txt Sun Apr 29 17:36:42 2018 -0400
b
@@ -1,7 +1,7 @@
 EntrezID Symbol GeneName Chr Length
+11302 Aatk apoptosis-associated tyrosine kinase 11 5743
 11287 Pzp pregnancy zone protein 6 4681
-11298 Aanat arylalkylamine N-acetyltransferase 11 1455
-11302 Aatk apoptosis-associated tyrosine kinase 11 5743
 11303 Abca1 ATP-binding cassette, sub-family A (ABC1), member 1 4 10260
 11304 Abca4 ATP-binding cassette, sub-family A (ABC1), member 4 3 7248
-11305 Abca2 ATP-binding cassette, sub-family A (ABC1), member 2 2 8061
\ No newline at end of file
+11305 Abca2 ATP-binding cassette, sub-family A (ABC1), member 2 2 8061
+11298 Aanat arylalkylamine N-acetyltransferase 11 1455
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-trend_Mut-WT.tsv
--- a/test-data/limma-trend_Mut-WT.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"logFC" "AveExpr" "t" "P.Value" "adj.P.Val" "B"
-"11304" 0.454006866294872 15.5347733798565 5.73783569533779 6.51380077190898e-06 3.90828046314539e-05 9.50615253604686
-"11287" 0.189211879829905 17.6560194022488 4.52571052740792 0.000138730279836171 0.000416190839508512 3.37443017677096
-"11298" -0.137906115120717 17.6760627799956 -3.28370943058883 0.00313397342215802 0.00626794684431604 -1.40612492217495
-"11303" -0.0565048762006306 17.887811868748 -1.22181580605621 0.233642901495507 0.35046435224326 -5.9848993971045
-"11305" -0.0597744562444475 18.1592280697576 -1.04836609029276 0.304913314924737 0.365895977909685 -6.17897539114193
-"11302" -0.0455235218343741 10.2501491061291 -0.373618198821297 0.711968519770644 0.711968519770644 -6.65188039014484
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-voom_Mut-WT.tsv
--- a/test-data/limma-voom_Mut-WT.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"GeneID" "logFC" "AveExpr" "t" "P.Value" "adj.P.Val" "B"
-"11304" 0.457332061341022 15.5254133001226 6.50459574633635 9.98720685007144e-07 5.99232411004286e-06 14.0741948485867
-"11287" 0.190749727701785 17.6546448244617 5.09535410066427 3.26518807653921e-05 9.79556422961764e-05 5.46773893802514
-"11298" -0.138014418336201 17.6747285193431 -3.3316848584234 0.00278753263633103 0.00557506527266206 -1.84301342041422
-"11303" -0.0558958943607024 17.886791401216 -1.30108531275582 0.205582481502277 0.254491025872966 -6.49241240578
-"11305" -0.0606991650996669 18.1585474109909 -1.28203791127301 0.212075854894138 0.254491025872966 -6.42090197700496
-"11302" -0.035023968220445 9.78883119065989 -0.236945963165271 0.814709535394086 0.814709535394086 -6.09497670655939
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-voom_Mut-WT_2fact.tsv
--- a/test-data/limma-voom_Mut-WT_2fact.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"GeneID" "logFC" "AveExpr" "t" "P.Value" "adj.P.Val" "B"
-"11287" 0.190521683937302 17.6546448244617 14.4538942795641 5.93483727833391e-09 3.56090236700035e-08 96.5864491179047
-"11298" -0.13987741320075 17.6747285193431 -7.7190169464239 5.41139783114199e-06 1.6234193493426e-05 22.3495333961476
-"11304" 0.459011254726059 15.5254133001226 5.6408319840902 0.00010884906092378 0.000217698121847559 8.76657226635715
-"11303" -0.0641599901594248 17.886791401216 -2.9900818153928 0.0112725655419901 0.0169088483129851 -2.70089035288878
-"11305" -0.0651479753016169 18.1585474109909 -2.28935282063866 0.0409794711446066 0.0491753653735279 -4.2713352660451
-"11302" -0.0358817134644287 9.78883119065989 -0.439436789626843 0.668154169055813 0.668154169055813 -5.75838315483931
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-voom_Mut-WT_2fact_anno.tsv
--- a/test-data/limma-voom_Mut-WT_2fact_anno.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"EntrezID" "Symbol" "GeneName" "Chr" "Length" "logFC" "AveExpr" "t" "P.Value" "adj.P.Val" "B"
-11287 "Pzp" "pregnancy zone protein" 6 4681 0.190521683937302 17.6546448244617 15.1152647604854 3.5605743515492e-09 2.13634461092952e-08 106.319792304784
-11298 "Aanat" "arylalkylamine N-acetyltransferase" 11 1455 -0.139877413200757 17.6747285193431 -8.49164953574146 2.03124167559177e-06 6.0937250267753e-06 28.5167891813121
-11304 "Abca4" "ATP-binding cassette, sub-family A (ABC1), member 4" 3 7248 0.459011254726057 15.5254133001226 5.60316118539212 0.000115556485215252 0.000231112970430504 8.56249981126828
-11303 "Abca1" "ATP-binding cassette, sub-family A (ABC1), member 1" 4 10260 -0.0641599901594141 17.886791401216 -2.92798536696768 0.0126512435710263 0.0189768653565394 -2.86386681629442
-11305 "Abca2" "ATP-binding cassette, sub-family A (ABC1), member 2" 2 8061 -0.0651479753016169 18.1585474109909 -2.03911882030581 0.0640903140266962 0.0769083768320354 -4.69855556350315
-11302 "Aatk" "apoptosis-associated tyrosine kinase" 11 5743 -0.0358817134644287 9.78883119065989 -0.43940476739557 0.668176726641997 0.668176726641997 -5.75841999159995
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-voom_Mut-WT_anno.tsv
--- a/test-data/limma-voom_Mut-WT_anno.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"EntrezID" "Symbol" "GeneName" "Chr" "Length" "logFC" "AveExpr" "t" "P.Value" "adj.P.Val" "B"
-11304 "Abca4" "ATP-binding cassette, sub-family A (ABC1), member 4" 3 7248 0.457332061341022 15.5254133001226 6.50459574633635 9.98720685007144e-07 5.99232411004286e-06 14.0741948485867
-11287 "Pzp" "pregnancy zone protein" 6 4681 0.190749727701785 17.6546448244617 5.09535410066427 3.26518807653921e-05 9.79556422961764e-05 5.46773893802514
-11298 "Aanat" "arylalkylamine N-acetyltransferase" 11 1455 -0.138014418336201 17.6747285193431 -3.3316848584234 0.00278753263633103 0.00557506527266206 -1.84301342041422
-11303 "Abca1" "ATP-binding cassette, sub-family A (ABC1), member 1" 4 10260 -0.0558958943607024 17.886791401216 -1.30108531275582 0.205582481502277 0.254491025872966 -6.49241240578
-11305 "Abca2" "ATP-binding cassette, sub-family A (ABC1), member 2" 2 8061 -0.0606991650996669 18.1585474109909 -1.28203791127301 0.212075854894138 0.254491025872966 -6.42090197700496
-11302 "Aatk" "apoptosis-associated tyrosine kinase" 11 5743 -0.035023968220445 9.78883119065989 -0.236945963165271 0.814709535394086 0.814709535394086 -6.09497670655939
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-voom_WT-Mut.tsv
--- a/test-data/limma-voom_WT-Mut.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"GeneID" "logFC" "AveExpr" "t" "P.Value" "adj.P.Val" "B"
-"11304" -0.457332061341022 15.5254133001226 -6.50459574633635 9.98720685007144e-07 5.99232411004286e-06 14.0741948485867
-"11287" -0.190749727701785 17.6546448244617 -5.09535410066427 3.26518807653921e-05 9.79556422961764e-05 5.46773893802514
-"11298" 0.138014418336201 17.6747285193431 3.3316848584234 0.00278753263633103 0.00557506527266206 -1.84301342041422
-"11303" 0.0558958943607024 17.886791401216 1.30108531275582 0.205582481502277 0.254491025872966 -6.49241240578
-"11305" 0.0606991650996669 18.1585474109909 1.28203791127301 0.212075854894138 0.254491025872966 -6.42090197700496
-"11302" 0.035023968220445 9.78883119065989 0.236945963165271 0.814709535394086 0.814709535394086 -6.09497670655939
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-voom_WT-Mut_2fact_anno.tsv
--- a/test-data/limma-voom_WT-Mut_2fact_anno.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"EntrezID" "Symbol" "GeneName" "Chr" "Length" "logFC" "AveExpr" "t" "P.Value" "adj.P.Val" "B"
-11287 "Pzp" "pregnancy zone protein" 6 4681 -0.190521683937302 17.6546448244617 -15.1152647604854 3.5605743515492e-09 2.13634461092952e-08 106.319792304784
-11298 "Aanat" "arylalkylamine N-acetyltransferase" 11 1455 0.139877413200757 17.6747285193431 8.49164953574146 2.03124167559177e-06 6.0937250267753e-06 28.5167891813121
-11304 "Abca4" "ATP-binding cassette, sub-family A (ABC1), member 4" 3 7248 -0.459011254726057 15.5254133001226 -5.60316118539212 0.000115556485215252 0.000231112970430504 8.56249981126828
-11303 "Abca1" "ATP-binding cassette, sub-family A (ABC1), member 1" 4 10260 0.0641599901594141 17.886791401216 2.92798536696768 0.0126512435710263 0.0189768653565394 -2.86386681629442
-11305 "Abca2" "ATP-binding cassette, sub-family A (ABC1), member 2" 2 8061 0.0651479753016169 18.1585474109909 2.03911882030581 0.0640903140266962 0.0769083768320354 -4.69855556350315
-11302 "Aatk" "apoptosis-associated tyrosine kinase" 11 5743 0.0358817134644287 9.78883119065989 0.43940476739557 0.668176726641997 0.668176726641997 -5.75841999159995
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-voom_normcounts.tsv
--- a/test-data/limma-voom_normcounts.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"GeneID" "Mut1" "Mut2" "Mut3" "WT1" "WT2" "WT3"
-"11287" 17.7719335903598 17.7105244453357 17.7658460810258 17.6076545522668 17.5079905058358 17.5639197719462
-"11298" 17.6506532355915 17.5520012519588 17.6144324004081 17.7727137985373 17.6986875511432 17.7598828784201
-"11302" 9.71615817858834 9.91760904198188 9.67886550454371 9.57719962386619 10.0175525101992 9.82560228478005
-"11303" 17.8774046020864 17.7865502833631 17.911613462219 17.9125899785601 17.8773613214092 17.955228759658
-"11304" 15.7545783969022 15.8519803740125 15.6561454280436 15.354518172169 15.2178020484983 15.3174553811097
-"11305" 18.0401341021246 18.1408682071359 18.204920085011 18.1808259688524 18.1818679259847 18.2026681768366
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/limma-voom_normcounts_anno.tsv
--- a/test-data/limma-voom_normcounts_anno.tsv Wed Jan 31 12:45:42 2018 -0500
+++ /dev/null Thu Jan 01 00:00:00 1970 +0000
b
@@ -1,7 +0,0 @@
-"EntrezID" "Symbol" "GeneName" "Chr" "Length" "Mut1" "Mut2" "Mut3" "WT1" "WT2" "WT3"
-11287 "Pzp" "pregnancy zone protein" 6 4681 17.7719335903598 17.7105244453357 17.7658460810258 17.6076545522668 17.5079905058358 17.5639197719462
-11298 "Aanat" "arylalkylamine N-acetyltransferase" 11 1455 17.6506532355915 17.5520012519588 17.6144324004081 17.7727137985373 17.6986875511432 17.7598828784201
-11302 "Aatk" "apoptosis-associated tyrosine kinase" 11 5743 9.71615817858834 9.91760904198188 9.67886550454371 9.57719962386619 10.0175525101992 9.82560228478005
-11303 "Abca1" "ATP-binding cassette, sub-family A (ABC1), member 1" 4 10260 17.8774046020864 17.7865502833631 17.911613462219 17.9125899785601 17.8773613214092 17.955228759658
-11304 "Abca4" "ATP-binding cassette, sub-family A (ABC1), member 4" 3 7248 15.7545783969022 15.8519803740125 15.6561454280436 15.354518172169 15.2178020484983 15.3174553811097
-11305 "Abca2" "ATP-binding cassette, sub-family A (ABC1), member 2" 2 8061 18.0401341021246 18.1408682071359 18.204920085011 18.1808259688524 18.1818679259847 18.2026681768366
b
diff -r 38aab66ae5cb -r a61a6e62e91f test-data/out_rscript.txt
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/out_rscript.txt Sun Apr 29 17:36:42 2018 -0400
[
b'@@ -0,0 +1,869 @@\n+# This tool takes in a matrix of feature counts as well as gene annotations and\n+# outputs a table of top expressions as well as various plots for differential\n+# expression analysis\n+#\n+# ARGS: htmlPath", "R", 1, "character"      -Path to html file linking to other outputs\n+#       outPath", "o", 1, "character"       -Path to folder to write all output to\n+#       filesPath", "j", 2, "character"     -JSON list object if multiple files input\n+#       matrixPath", "m", 2, "character"    -Path to count matrix\n+#       factFile", "f", 2, "character"      -Path to factor information file\n+#       factInput", "i", 2, "character"     -String containing factors if manually input\n+#       annoPath", "a", 2, "character"      -Path to input containing gene annotations\n+#       contrastData", "C", 1, "character"  -String containing contrasts of interest\n+#       cpmReq", "c", 2, "double"           -Float specifying cpm requirement\n+#       cntReq", "z", 2, "integer"          -Integer specifying minimum total count requirement\n+#       sampleReq", "s", 2, "integer"       -Integer specifying cpm requirement\n+#       normCounts", "x", 0, "logical"      -String specifying if normalised counts should be output\n+#       rdaOpt", "r", 0, "logical"          -String specifying if RData should be output\n+#       lfcReq", "l", 1, "double"           -Float specifying the log-fold-change requirement\n+#       pValReq", "p", 1, "double"          -Float specifying the p-value requirement\n+#       pAdjOpt", "d", 1, "character"       -String specifying the p-value adjustment method\n+#       normOpt", "n", 1, "character"       -String specifying type of normalisation used\n+#       robOpt", "b", 0, "logical"          -String specifying if robust options should be used\n+#       trend", "t", 1, "double"            -Float for prior.count if limma-trend is used instead of voom\n+#       weightOpt", "w", 0, "logical"       -String specifying if voomWithQualityWeights should be used\n+#\n+# OUT:\n+#       MDS Plot\n+#       Voom/SA plot\n+#       MD Plot\n+#       Expression Table\n+#       HTML file linking to the ouputs\n+# Optional:\n+#       Normalised counts Table\n+#       RData file\n+#\n+#\n+# Author: Shian Su - registertonysu@gmail.com - Jan 2014\n+# Modified by: Maria Doyle - Jun 2017, Jan 2018\n+\n+# Record starting time\n+timeStart <- as.character(Sys.time())\n+\n+# Load all required libraries\n+library(methods, quietly=TRUE, warn.conflicts=FALSE)\n+library(statmod, quietly=TRUE, warn.conflicts=FALSE)\n+library(splines, quietly=TRUE, warn.conflicts=FALSE)\n+library(edgeR, quietly=TRUE, warn.conflicts=FALSE)\n+library(limma, quietly=TRUE, warn.conflicts=FALSE)\n+library(scales, quietly=TRUE, warn.conflicts=FALSE)\n+library(getopt, quietly=TRUE, warn.conflicts=FALSE)\n+\n+if (packageVersion("limma") < "3.20.1") {\n+    stop("Please update \'limma\' to version >= 3.20.1 to run this tool")\n+}\n+\n+################################################################################\n+### Function Delcaration\n+################################################################################\n+# Function to sanitise contrast equations so there are no whitespaces\n+# surrounding the arithmetic operators, leading or trailing whitespace\n+sanitiseEquation <- function(equation) {\n+    equation <- gsub(" *[+] *", "+", equation)\n+    equation <- gsub(" *[-] *", "-", equation)\n+    equation <- gsub(" *[/] *", "/", equation)\n+    equation <- gsub(" *[*] *", "*", equation)\n+    equation <- gsub("^\\\\s+|\\\\s+$", "", equation)\n+    return(equation)\n+}\n+\n+# Function to sanitise group information\n+sanitiseGroups <- function(string) {\n+    string <- gsub(" *[,] *", ",", string)\n+    string <- gsub("^\\\\s+|\\\\s+$", "", string)\n+    return(string)\n+}\n+\n+# Function to change periods to whitespace in a string\n+unmake.names <- function(string) {\n+    string <- gsub(".", " ", string, fixed=TRUE)\n+    return(string)\n+}\n+\n+# Generate output folder and paths\n+makeOut <- function(filename) {\n+    return(past'..b'("<a href=\\"",\n+                  "http://www.bioconductor.org/packages/release/bioc/",\n+                  "vignettes/edgeR/inst/doc/edgeRUsersGuide.pdf",\n+                  "\\">", "edgeR User\'s Guide", "</a>")\n+\n+cit[1] <- paste("Please cite the following paper for this tool:")\n+\n+cit[2] <- paste("Liu R, Holik AZ, Su S, Jansz N, Chen K, Leong HS, Blewitt ME,",\n+                "Asselin-Labat ML, Smyth GK, Ritchie ME (2015). Why weight? ",\n+                "Modelling sample and observational level variability improves power ",\n+                "in RNA-seq analyses. Nucleic Acids Research, 43(15), e97.")\n+\n+cit[3] <- paste("Please cite the paper below for the limma software itself.",\n+                "Please also try to cite the appropriate methodology articles",\n+                "that describe the statistical methods implemented in limma,",\n+                "depending on which limma functions you are using. The",\n+                "methodology articles are listed in Section 2.1 of the",\n+                link[1],\n+                "Cite no. 3 only if sample weights were used.")\n+cit[4] <- paste("Smyth GK (2005). Limma: linear models for microarray data.",\n+                "In: \'Bioinformatics and Computational Biology Solutions using",\n+                "R and Bioconductor\'. R. Gentleman, V. Carey, S. doit,.",\n+                "Irizarry, W. Huber (eds), Springer, New York, pages 397-420.")\n+cit[5] <- paste("Please cite the first paper for the software itself and the",\n+                "other papers for the various original statistical methods",\n+                "implemented in edgeR.  See Section 1.2 in the", link[2],\n+                "for more detail.")\n+cit[6] <- paste("Robinson MD, McCarthy DJ and Smyth GK (2010). edgeR: a",\n+                "Bioconductor package for differential expression analysis",\n+                "of digital gene expression data. Bioinformatics 26, 139-140")\n+cit[7] <- paste("Robinson MD and Smyth GK (2007). Moderated statistical tests",\n+                "for assessing differences in tag abundance. Bioinformatics",\n+                "23, 2881-2887")\n+cit[8] <- paste("Robinson MD and Smyth GK (2008). Small-sample estimation of",\n+                "negative binomial dispersion, with applications to SAGE data.",\n+                "Biostatistics, 9, 321-332")\n+cit[9] <- paste("McCarthy DJ, Chen Y and Smyth GK (2012). Differential",\n+                "expression analysis of multifactor RNA-Seq experiments with",\n+                "respect to biological variation. Nucleic Acids Research 40,",\n+                "4288-4297")\n+cit[10] <- paste("Law CW, Chen Y, Shi W, and Smyth GK (2014). Voom:",\n+                "precision weights unlock linear model analysis tools for",\n+                "RNA-seq read counts. Genome Biology 15, R29.")\n+cit[11] <- paste("Ritchie ME, Diyagama D, Neilson J, van Laar R,",\n+                "Dobrovic A, Holloway A and Smyth GK (2006).",\n+                "Empirical array quality weights for microarray data.",\n+                "BMC Bioinformatics 7, Article 261.")\n+cata("<h3>Citations</h3>\\n")\n+cata(cit[1], "\\n")\n+cata("<br>\\n")\n+cata(cit[2], "\\n")\n+\n+cata("<h4>limma</h4>\\n")\n+cata(cit[3], "\\n")\n+cata("<ol>\\n")\n+ListItem(cit[4])\n+ListItem(cit[10])\n+ListItem(cit[11])\n+cata("</ol>\\n")\n+\n+cata("<h4>edgeR</h4>\\n")\n+cata(cit[5], "\\n")\n+cata("<ol>\\n")\n+ListItem(cit[6])\n+ListItem(cit[7])\n+ListItem(cit[8])\n+ListItem(cit[9])\n+cata("</ol>\\n")\n+\n+cata("<p>Please report problems or suggestions to: su.s@wehi.edu.au</p>\\n")\n+\n+for (i in 1:nrow(linkData)) {\n+    if (grepl("session_info", linkData$Link[i])) {\n+        HtmlLink(linkData$Link[i], linkData$Label[i])\n+    }\n+}\n+\n+cata("<table border=\\"0\\">\\n")\n+cata("<tr>\\n")\n+TableItem("Task started at:"); TableItem(timeStart)\n+cata("</tr>\\n")\n+cata("<tr>\\n")\n+TableItem("Task ended at:"); TableItem(timeEnd)\n+cata("</tr>\\n")\n+cata("<tr>\\n")\n+TableItem("Task run time:"); TableItem(timeTaken)\n+cata("<tr>\\n")\n+cata("</table>\\n")\n+\n+cata("</body>\\n")\n+cata("</html>")\n'