Files
galaxy/tools/maf/maf_to_interval.py
T

87 lines
2.8 KiB
Python

#!/usr/bin/env python
"""
Read a maf and output intervals for specified list of species.
"""
import os
import sys
from bx.align import maf
from galaxy.tools.util import maf_utilities
def __main__():
input_filename = sys.argv[1]
output_filename = sys.argv[2]
output_id = sys.argv[3]
# where to store files that become additional output
database_tmp_dir = sys.argv[4]
primary_spec = sys.argv[5]
species = sys.argv[6].split(",")
all_species = sys.argv[7].split(",")
partial = sys.argv[8]
keep_gaps = sys.argv[9]
out_files = {}
if "None" in species:
species = []
if primary_spec not in species:
species.append(primary_spec)
if primary_spec not in all_species:
all_species.append(primary_spec)
all_species.sort()
for spec in species:
if spec == primary_spec:
out_files[spec] = open(output_filename, "w+")
else:
out_files[spec] = open(
os.path.join(database_tmp_dir, "primary_%s_%s_visible_interval_%s" % (output_id, spec, spec)), "w+"
)
out_files[spec].write("#chrom\tstart\tend\tstrand\tscore\tname\t%s\n" % ("\t".join(all_species)))
num_species = len(all_species)
file_in = open(input_filename)
maf_reader = maf.Reader(file_in)
for i, m in enumerate(maf_reader):
for j, block in enumerate(maf_utilities.iter_blocks_split_by_species(m)):
if len(block.components) < num_species and partial == "partial_disallowed":
continue
sequences = {}
for c in block.components:
spec, chrom = maf_utilities.src_split(c.src)
if keep_gaps == "remove_gaps":
sequences[spec] = c.text.replace("-", "")
else:
sequences[spec] = c.text
sequences = "\t".join(sequences.get(_, "") for _ in all_species)
for spec in species:
c = block.get_component_by_src_start(spec)
if c is not None:
spec2, chrom = maf_utilities.src_split(c.src)
assert spec2 == spec, Exception(
"Species name inconsistancy found in component: %s != %s" % (spec, spec2)
)
out_files[spec].write(
"%s\t%s\t%s\t%s\t%s\t%s\t%s\n"
% (
chrom,
c.forward_strand_start,
c.forward_strand_end,
c.strand,
m.score,
"%s_%s_%s" % (spec, i, j),
sequences,
)
)
file_in.close()
for file_out in out_files.values():
file_out.close()
if __name__ == "__main__":
__main__()