Mercurial > repos > dereeper > pangenome_explorer
comparison Snakemake_files/Snakefile_wget_roary_heatmap_upset_COG @ 3:e42d30da7a74 draft
Uploaded
author | dereeper |
---|---|
date | Thu, 30 May 2024 11:52:25 +0000 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
2:97e4e3e818b6 | 3:e42d30da7a74 |
---|---|
1 import glob | |
2 import os | |
3 import shutil | |
4 import re | |
5 import subprocess | |
6 | |
7 import yaml | |
8 configfile: "config.yaml" | |
9 | |
10 | |
11 SAMPLES = [] | |
12 with open("config.yaml", "r") as yaml_file: | |
13 genome_data = yaml.safe_load(yaml_file) | |
14 if "ids" in genome_data.keys(): | |
15 for id in genome_data["ids"]: | |
16 SAMPLES.append(id) | |
17 if "input_genbanks" in genome_data.keys(): | |
18 for gb_path in genome_data["input_genbanks"]: | |
19 | |
20 cmd = "grep 'ACCESSION' "+gb_path | |
21 returned_value = subprocess.getoutput(cmd) | |
22 words = returned_value.split() | |
23 name = words[1] | |
24 #words = returned_value.split() | |
25 #SAMPLES.append(words[1]) | |
26 words = gb_path.split("/") | |
27 genbank_file_name = words[-1] | |
28 x = genbank_file_name.replace(".", "_") | |
29 cmd = "sed -i 's/ACCESSION "+name+"/ACCESSION "+x+"/g' "+gb_path | |
30 subprocess.getoutput(cmd) | |
31 | |
32 cmd_locus = "grep 'LOCUS' "+gb_path | |
33 returned_value = subprocess.getoutput(cmd_locus) | |
34 words_locus = returned_value.split() | |
35 name_locus = words_locus[1] | |
36 cmd = "sed -i 's/LOCUS "+name_locus+"/LOCUS "+x+"/' "+gb_path | |
37 subprocess.getoutput(cmd) | |
38 | |
39 cmd = "grep 'ACCESSION' "+gb_path | |
40 returned_value = subprocess.getoutput(cmd) | |
41 words = returned_value.split() | |
42 SAMPLES.append(words[1]) | |
43 | |
44 | |
45 rule final: | |
46 input: | |
47 "outputs/pav_matrix.csv", | |
48 "outputs/GCskew.txt", | |
49 "outputs/pav_matrix.tsv", | |
50 "outputs/heatmap.svg.gz", | |
51 "outputs/cog_output.txt", | |
52 "outputs/fastani.out", | |
53 "outputs/rarefaction_curves.txt" | |
54 | |
55 | |
56 rule ncbi_datasets: | |
57 input: | |
58 "config.yaml" | |
59 output: | |
60 expand("outputs/genomes/{sample}.fasta", sample=SAMPLES), | |
61 expand("outputs/genomes/{sample}.gb", sample=SAMPLES), | |
62 genomes="outputs/genomes/genomes.txt" | |
63 shell: | |
64 """ | |
65 perl $PANEX_PATH/Perl/get_data.pl {input} outputs/genomes | |
66 """ | |
67 | |
68 rule gcskew: | |
69 input: | |
70 "outputs/genomes/{sample}.fasta" | |
71 output: | |
72 "outputs/genomes/{sample}.fasta.gcskew.txt" | |
73 shell: | |
74 """ | |
75 python3 $PANEX_PATH/SkewIT/src/gcskew.py -i {input} -o {input}.gcskew.txt -k 1000 -w 1000 | |
76 """ | |
77 | |
78 rule concat_gcskew: | |
79 input: | |
80 expand("outputs/genomes/{sample}.fasta.gcskew.txt", sample=SAMPLES) | |
81 output: | |
82 out2="outputs/GCskew.txt" | |
83 shell: | |
84 """ | |
85 cat {input} >>{output.out2} | |
86 """ | |
87 | |
88 rule genbank2gff3: | |
89 input: | |
90 "outputs/genomes/{sample}.gb" | |
91 output: | |
92 gff1="outputs/genomes/{sample}.gb.gff", | |
93 gff2="outputs/genomes/{sample}.gb.rmdup.gff" | |
94 shell: | |
95 """ | |
96 perl $PANEX_PATH/Perl/bp_genbank2gff3.pl -o outputs/genomes {input} | |
97 perl $PANEX_PATH/Perl/remove_duplicates_in_gff.pl {output.gff1} {output.gff2} | |
98 """ | |
99 | |
100 rule roary: | |
101 input: | |
102 expand("outputs/genomes/{sample}.gb.rmdup.gff", sample=SAMPLES) | |
103 params: | |
104 identity=80 | |
105 output: | |
106 pav="outputs/pav_matrix.csv", | |
107 newick="outputs/accessory_binary_genes.fa.newick" | |
108 shell: | |
109 """ | |
110 roary -f outputs/roary_outdir -i {params.identity} outputs/genomes/*gb.rmdup.gff | |
111 cp -rf outputs/roary_outdir/gene_presence_absence.csv {output.pav} | |
112 cp -rf outputs/roary_outdir/accessory_binary_genes.fa.newick {output.newick} | |
113 """ | |
114 | |
115 rule convert_matrix: | |
116 input: | |
117 pav="outputs/pav_matrix.csv" | |
118 output: | |
119 "outputs/pav_matrix.tsv" | |
120 shell: | |
121 """ | |
122 perl $PANEX_PATH/Perl/ConvertRoaryMatrix.pl outputs/genomes {input} {output} outputs/genomes/strains.txt | |
123 """ | |
124 | |
125 rule cog: | |
126 input: | |
127 pav="outputs/pav_matrix.tsv" | |
128 output: | |
129 cog="outputs/cog_output.txt", | |
130 cogstat="outputs/cog_stats.txt", | |
131 cogstat2="outputs/cog_stats2.txt", | |
132 cogofclusters="outputs/cog_of_clusters.txt" | |
133 shell: | |
134 """ | |
135 perl $PANEX_PATH/Perl/GetCogOfCluster.pl {input} outputs/genomes {output.cog} {output.cogstat} {output.cogstat2} {output.cogofclusters} outputs/genomes/strains.txt | |
136 """ | |
137 | |
138 rule heatmap_upset: | |
139 input: | |
140 pav="outputs/pav_matrix.tsv" | |
141 output: | |
142 heatmap="outputs/heatmap.svg.gz", | |
143 html="outputs/heatmap.svg.heatmap_plotly.html", | |
144 upsetr="outputs/upsetr.svg", | |
145 binpav="outputs/heatmap.svg.pangenes_01matrix.txt" | |
146 shell: | |
147 """ | |
148 perl $PANEX_PATH/Perl/GenerateHeatmapFromPAV.pl {input.pav} outputs/heatmap.svg | |
149 mv outputs/heatmap.svg.upsetr.svg {output.upsetr} | |
150 """ | |
151 | |
152 rule micropan: | |
153 input: | |
154 binpav="outputs/heatmap.svg.pangenes_01matrix.txt" | |
155 output: | |
156 txt="outputs/rarefaction_curves.txt", | |
157 pdf="outputs/rarefaction_curves.pdf", | |
158 svg="outputs/rarefaction_curves.svg", | |
159 heaps="outputs/heaps.tsv" | |
160 shell: | |
161 """ | |
162 Rscript $PANEX_PATH/R/micropan_rarefaction.R -f {input.binpav} -p {output.pdf} -a {output.heaps} -o {output.txt} | |
163 pdf2svg {output.pdf} {output.svg} | |
164 """ | |
165 | |
166 rule fastani: | |
167 input: | |
168 "outputs/genomes/genomes.txt" | |
169 output: | |
170 out="outputs/fastani.out", | |
171 matrix="outputs/fastani.out.matrix", | |
172 completematrix="outputs/fastani.out.matrix.complete", | |
173 pdf="outputs/fastani.out.pdf", | |
174 svg="outputs/fastani.out.svg" | |
175 shell: | |
176 """ | |
177 fastANI --rl {input} --ql {input} -o {output.out} -t 4 --matrix | |
178 perl $PANEX_PATH/Perl/convertANI.pl {output.matrix} outputs/genomes/genomes2.txt >{output.completematrix} | |
179 Rscript $PANEX_PATH/R/heatmap_ani.R -f {output.completematrix} -o {output.pdf} | |
180 pdf2svg {output.pdf} {output.svg} | |
181 """ | |
182 |