https://github.com/fjruizruano/ngs-protocols
Tip revision: 39a091d1fa569a7fc717ac73c4b3de07f0a1204d authored by fjruizruano on 03 August 2023, 11:48:27 UTC
adding gfa2fas.py and extract_gfa.py
adding gfa2fas.py and extract_gfa.py
Tip revision: 39a091d
coverage_window.py
#!/usr/bin/python
import sys
from subprocess import call
print "Usage: coverage_window.py BedGraphFile FastaFile WindowSize"
try:
bg = sys.argv[1]
except:
bg = raw_input("Introduce BedGraph file: ")
try:
ref = sys.argv[2]
except:
ref = raw_input("Introduce Fasta file: ")
try:
window = sys.argv[3]
window = int(window)
except:
window = raw_input("Introduce window size (integer): ")
window = int(window)
data = open(bg).readlines()
header = data[0]
header = header.split()
samples = header[3:]
try:
fai = open(ref+".fai").readlines()
except:
call("samtools faidx %s" % ref, shell=True)
fai = open(ref+".fai").readlines()
lis = [0] * len(samples)
results = {}
for line in fai:
line = line.split()
name = line[0]
length = int(line[1])
results[name] = {}
for n in range(0,(length/window)+1):
results[name][n] = lis
for line in data[1:]:
line = line.split()
name = line[0]
start = int(line[1])
end = int(line[2])
count = [float(i) for i in line[3:]]
for n in range(0,end-start):
number = start+n
count_old = results[name][number/window]
count_new = [x+y for x, y in zip(count_old, count)]
results[name][number/window] = count_new
w = open(bg+".counts"+str(window), "w")
w.write("Contig\tWindow\t%s\n" % "\t".join(samples))
for contig in results:
for win in results[contig]:
info = results[contig][win]
info = [str(int(i)) for i in info]
w.write("%s\t%s\t%s\n" % (contig,str(win),"\t".join(info)))