comparison tools/protein_analysis/promoter2.py @ 20:a19b3ded8f33 draft

v0.2.11 Job splitting fast-fail; RXLR tools supports HMMER2 from BioConda; Capture more version information; misc internal changes

19:f3ecd80850e2

tab separated table.

Additionally, in order to take advantage of multiple cores the input FASTA
file is broken into chunks and multiple copies of promoter run at once.
This can be used in combination with the job-splitting available in Galaxy.

Note that rewriting the FASTA input file allows us to avoid a bug in
promoter 2 with long descriptions in the FASTA header line (over 200
characters) which produces stray fragements of the description in the
output file, making parsing non-trivial.

TODO - Automatically extract the sequence containing a promoter prediction?
"""
import sys
import os
import commands
import tempfile
from seq_analysis_utils import split_fasta, run_jobs, thread_count

FASTA_CHUNK = 500

if "-v" in sys.argv or "--version" in sys.argv:
    sys.exit(os.system("promoter -V"))

tmp_dir = tempfile.mkdtemp()


def get_path_and_binary():
    platform = commands.getoutput("uname")  # e.g. Linux
    shell_script = commands.getoutput("which promoter")
    if not os.path.isfile(shell_script):
        sys.exit("ERROR: Missing promoter executable shell script")
    path = None

def make_tabular(raw_handle, out_handle):
    """Parse text output into tabular, return query count."""
    identifier = None
    queries = 0
    for line in raw_handle:
        # print repr(line)
        if not line.strip() or line == "Promoter prediction:\n":
            pass
        elif line[0] != " ":
            identifier = line.strip().replace("\t", " ").split(None, 1)[0]
            queries += 1

            assert identifier
        else:
            try:
                position, score, likelihood = line.strip().split(None, 2)
            except ValueError:
                print "WARNING: Problem with line: %r" % line
                continue
                # sys.exit("ERROR: Problem with line: %r" % line)
            if likelihood not in ["ignored",
                                  "Marginal prediction",
                                  "Medium likely prediction",
                                  "Highly likely prediction"]:
                sys.exit("ERROR: Problem with line: %r" % line)
            out_handle.write("%s\t%s\t%s\t%s\n" % (identifier, position, score, likelihood))
    return queries

working_dir, bin = get_path_and_binary()

if not os.path.isfile(fasta_file):
    sys.exit("ERROR: Missing input FASTA file %r" % fasta_file)

    try:
        os.rmdir(tmp_dir)
    except Exception:
        pass

if len(jobs) > 1 and num_threads > 1:
    # A small "info" message for Galaxy to show the user.
    print "Using %i threads for %i tasks" % (min(num_threads, len(jobs)), len(jobs))
cur_dir = os.path.abspath(os.curdir)
os.chdir(working_dir)
results = run_jobs(jobs, num_threads)
os.chdir(cur_dir)
for fasta, temp, cmd in zip(fasta_files, temp_files, jobs):

        sys.exit("No output from promoter2")
    queries += count
out_handle.close()

clean_up(fasta_files + temp_files)
print "Results for %i queries" % queries

20:a19b3ded8f33

tab separated table.

Additionally, in order to take advantage of multiple cores the input FASTA
file is broken into chunks and multiple copies of promoter run at once.
This can be used in combination with the job-splitting available in Galaxy.
Note that rewriting the FASTA input file allows us to avoid a bug in
promoter 2 with long descriptions in the FASTA header line (over 200
characters) which produces stray fragements of the description in the
output file, making parsing non-trivial.

TODO - Automatically extract the sequence containing a promoter prediction?
"""

from __future__ import print_function

import commands
import os
import sys
import tempfile

from seq_analysis_utils import run_jobs, split_fasta, thread_count

FASTA_CHUNK = 500

if "-v" in sys.argv or "--version" in sys.argv:
    sys.exit(os.system("promoter -V"))

tmp_dir = tempfile.mkdtemp()


def get_path_and_binary():
    """Determine path and binary names for promoter tool."""
    platform = commands.getoutput("uname")  # e.g. Linux
    shell_script = commands.getoutput("which promoter")
    if not os.path.isfile(shell_script):
        sys.exit("ERROR: Missing promoter executable shell script")
    path = None

def make_tabular(raw_handle, out_handle):
    """Parse text output into tabular, return query count."""
    identifier = None
    queries = 0
    for line in raw_handle:
        # print(repr(line))
        if not line.strip() or line == "Promoter prediction:\n":
            pass
        elif line[0] != " ":
            identifier = line.strip().replace("\t", " ").split(None, 1)[0]
            queries += 1

            assert identifier
        else:
            try:
                position, score, likelihood = line.strip().split(None, 2)
            except ValueError:
                print("WARNING: Problem with line: %r" % line)
                continue
                # sys.exit("ERROR: Problem with line: %r" % line)
            if likelihood not in ["ignored",
                                  "Marginal prediction",
                                  "Medium likely prediction",
                                  "Highly likely prediction"]:
                sys.exit("ERROR: Problem with line: %r" % line)
            out_handle.write("%s\t%s\t%s\t%s\n" % (identifier, position, score, likelihood))
    return queries


working_dir, bin = get_path_and_binary()

if not os.path.isfile(fasta_file):
    sys.exit("ERROR: Missing input FASTA file %r" % fasta_file)

    try:
        os.rmdir(tmp_dir)
    except Exception:
        pass


if len(jobs) > 1 and num_threads > 1:
    # A small "info" message for Galaxy to show the user.
    print("Using %i threads for %i tasks" % (min(num_threads, len(jobs)), len(jobs)))
cur_dir = os.path.abspath(os.curdir)
os.chdir(working_dir)
results = run_jobs(jobs, num_threads)
os.chdir(cur_dir)
for fasta, temp, cmd in zip(fasta_files, temp_files, jobs):

        sys.exit("No output from promoter2")
    queries += count
out_handle.close()

clean_up(fasta_files + temp_files)
print("Results for %i queries" % queries)