Files
galaxy/tools/maf/maf_to_bed.py
T

119 lines
3.5 KiB
Python

#!/usr/bin/env python
"""
Read a maf and output intervals for specified list of species.
"""
from __future__ import print_function
import os
import sys
from collections import OrderedDict
from bx.align import maf
from bx.align.core import src_split
def __main__():
input_filename = sys.argv[1]
output_filename = sys.argv[2]
# where to store files that become additional output
database_tmp_dir = sys.argv[5]
species = sys.argv[3].split(",")
partial = sys.argv[4]
output_id = sys.argv[6]
out_files = {}
primary_spec = None
if "None" in species:
species = OrderedDict()
try:
for m in maf.Reader(open(input_filename)):
for c in m.components:
spec, chrom = src_split(c.src)
if not spec or not chrom:
spec = c.src
species[spec] = None
species = species.keys()
except Exception:
print("Invalid MAF file specified", file=sys.stderr)
return
if "?" in species:
print("Invalid dbkey specified", file=sys.stderr)
return
for i, spec in enumerate(species):
if i == 0:
out_files[spec] = open(output_filename, "w")
primary_spec = spec
else:
out_files[spec] = open(
os.path.join(database_tmp_dir, "primary_%s_%s_visible_bed_%s" % (output_id, spec, spec)), "w+"
)
num_species = len(species)
print("Restricted to species:", ",".join(species))
with open(input_filename) as file_in:
maf_reader = maf.Reader(file_in)
block_num = -1
for m in maf_reader:
block_num += 1
if "None" not in species:
m = m.limit_to_species(species)
components = m.components
if len(components) < num_species and partial == "partial_disallowed":
continue
for c in components:
spec, chrom = src_split(c.src)
if not spec or not chrom:
spec = chrom = c.src
if spec not in out_files.keys():
out_files[spec] = open(
os.path.join(database_tmp_dir, "primary_%s_%s_visible_bed_%s" % (output_id, spec, spec)), "wb+"
)
if c.strand == "-":
out_files[spec].write(
chrom
+ "\t"
+ str(c.src_size - c.end)
+ "\t"
+ str(c.src_size - c.start)
+ "\t"
+ spec
+ "_"
+ str(block_num)
+ "\t"
+ "0\t"
+ c.strand
+ "\n"
)
else:
out_files[spec].write(
chrom
+ "\t"
+ str(c.start)
+ "\t"
+ str(c.end)
+ "\t"
+ spec
+ "_"
+ str(block_num)
+ "\t"
+ "0\t"
+ c.strand
+ "\n"
)
for file_out in out_files.keys():
out_files[file_out].close()
print("#FILE1_DBKEY\t%s" % (primary_spec))
if __name__ == "__main__":
__main__()