Updated GO enrichment analysis to save plot.

Snakefile

@@ -17,7 +17,7 @@ include: "rules/ComparisonRegions.smk"
 include: "rules/Circos.smk"
 include: "rules/SampleSizes.smk"
 include: "rules/UseControl.smk"
-include: "rules/ComparisonSettings.smk"
+#include: "rules/ComparisonSettings.smk"
 
 #Directory to save all files.
 workdir: config["workdir"]
@@ -43,10 +43,10 @@ onerror:
 rule all:
     """Define desired output from pipeline."""
     input:
-        "GO/comparison.txt",
-        "Settings/Report.txt",
-        #"Samplesize/Report.txt",
         #"Reports/Tools.txt",
+        "GO/comparison.txt",
+        #"Settings/Report.txt",
+        #"Samplesize/Report.txt",aomt
         #"Samplesize/Report.txt"
 
 rule help:

config.yaml

@@ -21,7 +21,7 @@ biomart_genes: input_files/biomart_human_genes.tsv
 ID_to_GO: input_files/ID_to_GO.txt
 
 #Settings GISTIC2.0
-gistic_precision: "99"
+gistic_precision: "75"
 settings_gistic: "-ta 0.1
                   -td 0.1
                   -qvt 0.25
@@ -29,7 +29,7 @@ settings_gistic: "-ta 0.1
                   -cap 1.5
                   -rx 1
                   -genegistic 1
-                  -conf 0.99"
+                  -conf 0.75"
 
 comparison_settings: ["-ta 0.1 -td 0.1 -qvt 0.25 -brlen 0.7 -cap 1.5 -rx 1 -genegistic 1 -conf 0.99",
                       "-ta 0.1 -td 0.1 -qvt 0.25 -brlen 0.7 -cap 1.5 -rx 1 -genegistic 1 -conf 0.75",

rules/ComparisonSettings.smk

@@ -7,12 +7,12 @@ rule gistic_settings:
         seg="Input/Segments_tumor.txt"
     output:
         expand("Settings/GISTIC_{setting_nr}/all_lesions.conf_" + config["gistic_precision"] + ".txt", setting_nr=range(len(config["comparison_settings"]))),
-        expand("Settings/GISTIC_{setting_nr}/regions_track.conf_" + config["gistic_precision"] + ".bed", setting_nr=range(len(config["comparison_settings"])))
+        "Settings/GISTIC_{setting_nr}/regions_track.conf_" + config["gistic_precision"] + ".bed"
     params:
         cnv="",
         ref=config["reference"],
         ref_file="",
-        extra=expand("{setting}", setting=config["comparison_settings"])
+        extra=(config["comparison_settings"])[int("{setting_nr}")]
     wrapper:
         "file:" + workflow.basedir + "/wrappers/GISTIC2"
 

rules/GenePrioritization.smk

@@ -3,7 +3,8 @@ rule go_analysis:
     input:
         gene_list="Reports/Genes_{tool}.txt"
     output:
-        go="GO/{tool}.txt"
+        table="GO/Enriched_GOs_{tool}.txt",
+        plot="GO/Enriched_GOs_{tool}.jpg"
     params:
         organism="hsapiens",    #default is human
         ontology="BP",           #MF, BP, CC or all
@@ -16,7 +17,7 @@ rule go_analysis:
 rule compare_enriched_terms:
     """Compare the top 50 GO terms detected by RUBIC and GISTIC2.0"""
     input:
-        go=expand("GO/{tool}.txt", tool=["GISTIC2", "RUBIC"])
+        go=expand("GO/Enriched_GOs_{tool}.txt", tool=["GISTIC2", "RUBIC"])
     output:
         "GO/comparison.txt"
     run:

scripts/ReportTools.py

@@ -104,12 +104,13 @@ def venn3_overlap(gene_lists, known_genes, out_venn):
 def plot_histogram_sizes(parsed_results, plot_file):
     """Make boxplot with swarmplot on top showing sizes of recurrent regions detected by GISTIC2.0 and RUBIC."""
     size_list, size_list_known = parse_sizes_for_plot(parsed_results)
+    print(size_list_known)
     sns.set_style("white")
     f, ax = plt.subplots(figsize=(7, 7))
     ax.set(yscale="log")
     g = sns.boxplot(x="Tool", y="Size", hue="CNV type", data=size_list, whis=np.inf, palette={"Amplification": "#71AEC0", "Deletion": "#B55D60"})
     g = sns.swarmplot(x="Tool", y="Size", dodge=True, hue="CNV type", data=size_list, size=5, palette={"Amplification": "#527F8C", "Deletion": "#844446"})
-    g = sns.swarmplot(x="Tool", y="Size", dodge=True, hue="CNV type", data=size_list_known, size=5, palette={"Amplification": "black", "Deletion": "black"})
+    #g = sns.swarmplot(x="Tool", y="Size", dodge=True, hue="CNV type", data=size_list_known, size=5, palette={"Amplification": "black", "Deletion": "black"})
     sns.despine()
     g.set_xlabel("Tool",fontsize=16)
     g.set_ylabel("Log size of recurrent regions",fontsize=16)

wrappers/topgo/wrapper.R

@@ -48,5 +48,7 @@ print(allRes)
 write.table(allRes, file=snakemake@output[[1]], quote=FALSE, sep="\t", eol = "\n")
 
 #Visualize GO term relationships
+jpeg(snakemake@output[[1]])
 showSigOfNodes(GOdata, score(resultFisher), firstSigNodes = 5, useInfo ='all')
 printGraph(GOdata, resultFisher, firstSigNodes = 5, fn.prefix = "tGO", useInfo = "all", pdfSW = TRUE)
+dev.off()