find_meme_letter_probabilities.py

   1 #!/usr/bin/python 3
   2 import argparse
   3 import re
   4 import os
   5
   6 def main(cmdline=None):
   7     #filename = '/woldlab/castor/home/kfisher/ChIA/AliNodes/DNase4/Wellington/PP_DE_analysis/Consolidated_Wellington_promoters_3xUp_10FPKM.cyteConnected_myocyteFootprints.20rad.MEME2_maxw25/meme.txt'
   8     #output_dir = '.'
   9
  10     parser = make_parser()
  11     args = parser.parse_args(cmdline)
  12
  13     if not os.path.isdir(args.output):
  14         os.mkdir(args.output)
  15
  16     for filename in args.filenames:
  17         with open(filename, 'rt') as instream:
  18             version = find_meme_version(instream)
  19
  20             motif_counter = 0
  21             for motif in find_meme_letter_probabilities(instream):
  22                 motif_counter += 1
  23                 output_name = 'motif_%d.mot' % (motif_counter,)
  24                 output_filename = os.path.join(args.output, output_name)
  25                 with open(output_filename, 'wt') as outstream:
  26                     save_meme_letter_probabilities(
  27                         outstream, version, motif_counter, motif)
  28
  29
  30 def make_parser():
  31     parser = argparse.ArgumentParser('Search through a MEME report for letter-probability matrices')
  32     parser.add_argument('-o', '--output', default='.',
  33                         help="output directory")
  34     parser.add_argument('filenames', nargs='*',
  35                         help="name of MEME output files to search through")
  36     return parser
  37
  38 def find_meme_version(instream, limit=20):
  39     """Read through the first few meme header lines looking for the MEME version
  40     """
  41     version_pattern = re.compile("MEME version (?P<ver>[\d.]+)")
  42
  43     for i, line in enumerate(instream):
  44         if i > limit:
  45             return None
  46         else:
  47             match = version_pattern.match(line)
  48             if match:
  49                 return match.group('ver')
  50
  51
  52 def find_meme_letter_probabilities(instream):
  53     """Look through a meme output report and find the letter-probability matrix
  54
  55     Arguments:
  56       instream (input stream): file-like object to read from
  57
  58     Yields:
  59       list of strings for a MEME letter-probabilities matrix
  60     """
  61     header = "letter-probability matrix: alength=(?P<letters>[ \d]+) w=(?P<bases>[ \d]+)"
  62     letter_pattern = re.compile(header)
  63
  64     motif = []
  65     rows = 0
  66     for line in instream:
  67         line = line.rstrip()
  68         match = letter_pattern.match(line)
  69         if match:
  70             if motif:
  71                 yield motif
  72                 motif = []
  73             rows = int(match.group('bases'))
  74             motif.append(line)
  75         elif rows:
  76             motif.append(line)
  77             rows -= 1
  78     yield motif
  79
  80
  81 def save_meme_letter_probabilities(outstream, version, motif_index, motif):
  82     """Write MEME mot formatted contents to outstream
  83
  84     Arguments:
  85       outstream (output stream): file-like object to write to
  86       version (str): motif version number
  87       motif_index (int): what motif number we're on
  88       motif (list of strings): motif letter-probabilities from meme
  89     """
  90
  91     report = """MEME version {version}
  92
  93 ALPHABET= ACGT
  94
  95 strands: + -
  96
  97 Background etter frequencies (from uniform background):
  98 A 0.25000 C 0.25000 G 0.25000 T 0.25000
  99
 100 MOTIF {count}
 101
 102 {matrix}
 103 """.format(version=version,
 104            count=motif_index,
 105            matrix=os.linesep.join(motif))
 106     outstream.write(report)
 107
 108
 109 if __name__ == "__main__":
 110     main()
 111