find_meme_letter_probabilities.py

   1 #!/usr/bin/python 3
   2 """Program to search through MEME report to find letter-probabilties matrix
   3 """
   4
   5 # Copyright (2015) Diane Trout & California Institute of Technology
   6
   7 # This program is free software; you can redistribute it and/or modify
   8 # it under the terms of the GNU General Public License as published by
   9 # the Free Software Foundation; either version 2 of the License, or
  10 # (at your option) any later version.
  11 #
  12 # This program is distributed in the hope that it will be useful,
  13 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  14 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15 # GNU General Public License for more details.
  16 #
  17 # You should have received a copy of the GNU General Public License along
  18 # with this program; if not, write to the Free Software Foundation, Inc.,
  19 # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
  20 import argparse
  21 import re
  22 import os
  23
  24 def main(cmdline=None):
  25     parser = make_parser()
  26     args = parser.parse_args(cmdline)
  27
  28     if not os.path.isdir(args.output):
  29         os.mkdir(args.output)
  30
  31     for filename in args.filenames:
  32         with open(filename, 'rt') as instream:
  33             version = find_meme_version(instream)
  34
  35             motif_counter = 0
  36             for motif in find_meme_letter_probabilities(instream):
  37                 motif_counter += 1
  38                 output_name = 'motif_%d.mot' % (motif_counter,)
  39                 output_filename = os.path.join(args.output, output_name)
  40                 with open(output_filename, 'wt') as outstream:
  41                     save_meme_letter_probabilities(
  42                         outstream, version, motif_counter, motif)
  43
  44
  45 def make_parser():
  46     parser = argparse.ArgumentParser('Search through a MEME report for letter-probability matrices')
  47     parser.add_argument('-o', '--output', default='.',
  48                         help="output directory")
  49     parser.add_argument('filenames', nargs='*',
  50                         help="name of MEME output files to search through")
  51     return parser
  52
  53 def find_meme_version(instream, limit=20):
  54     """Read through the first few meme header lines looking for the MEME version
  55     """
  56     version_pattern = re.compile("MEME version (?P<ver>[\d.]+)")
  57
  58     for i, line in enumerate(instream):
  59         if i > limit:
  60             return None
  61         else:
  62             match = version_pattern.match(line)
  63             if match:
  64                 return match.group('ver')
  65
  66
  67 def find_meme_letter_probabilities(instream):
  68     """Look through a meme output report and find the letter-probability matrix
  69
  70     Arguments:
  71       instream (input stream): file-like object to read from
  72
  73     Yields:
  74       list of strings for a MEME letter-probabilities matrix
  75     """
  76     header = "letter-probability matrix: alength=(?P<letters>[ \d]+) w=(?P<bases>[ \d]+)"
  77     letter_pattern = re.compile(header)
  78
  79     motif = []
  80     rows = 0
  81     for line in instream:
  82         line = line.rstrip()
  83         match = letter_pattern.match(line)
  84         if match:
  85             if motif:
  86                 yield motif
  87                 motif = []
  88             rows = int(match.group('bases'))
  89             motif.append(line)
  90         elif rows:
  91             motif.append(line)
  92             rows -= 1
  93     yield motif
  94
  95
  96 def save_meme_letter_probabilities(outstream, version, motif_index, motif):
  97     """Write MEME mot formatted contents to outstream
  98
  99     Arguments:
 100       outstream (output stream): file-like object to write to
 101       version (str): motif version number
 102       motif_index (int): what motif number we're on
 103       motif (list of strings): motif letter-probabilities from meme
 104     """
 105
 106     report = """MEME version {version}
 107
 108 ALPHABET= ACGT
 109
 110 strands: + -
 111
 112 Background etter frequencies (from uniform background):
 113 A 0.25000 C 0.25000 G 0.25000 T 0.25000
 114
 115 MOTIF {count}
 116
 117 {matrix}
 118 """.format(version=version,
 119            count=motif_index,
 120            matrix=os.linesep.join(motif))
 121     outstream.write(report)
 122
 123
 124 if __name__ == "__main__":
 125     main()
 126