allpy: c1091715f8a3 test/usecase2.py

allpy

view test/usecase2.py @ 1088:c1091715f8a3

Added allpy.util.Silence to faciliate running most code with as little noise as possible This class is intended for use around noisy chunks of code using `with`, e.g.:: with Silence(): structure = PDBParser().load_pdb(file) This will produce no output to stderr if no exceptions are raised in the code. Otherwise, all output will be printed as if the code was run without this wrapper. See docstrings for more documentation.

author	Daniil Alexeyevsky <dendik@kodomo.fbb.msu.ru>
date	Fri, 01 Jun 2012 13:43:02 +0400
parents	cc1959669928
children	08d892230e8c

line source

1 #!/usr/bin/python

2 import sys

3 from collections import deque

4 from allpy import dna

5 from allpy.processors import Needle, Left

6 from allpy.fileio import FastaFile

7 from allpy.util import open

9 width = 15

10 threshold = 10

11 if __name__ == "__main__":

12 infile = open(sys.argv[1])

13 outfile = open(sys.argv[2])

15 def has_identity(column):

16 as_list = column.values()

17 return len(column) == 2 and as_list[0] == as_list[1]

19 def is_good_window(window):

20 sum_id = sum(int(has_identity(column)) for column in window)

21 return len(window) == width and sum_id >= threshold

23 def find_runs(alignment):

24 window = deque([], width)

25 blocks = []

26 in_block = False

27 for column in alignment.columns:

28 window.append(column)

29 in_block, was_in_block = is_good_window(window), in_block

30 if in_block and not was_in_block:

31 block = dna.Block.from_alignment(alignment, columns=list(window))

32 blocks.append(block)

33 elif in_block:

34 block.columns.append(column)

35 return blocks

37 def blocks_markup(alignment, blocks):

38 for column in alignment.columns:

39 column.in_block = "-"

40 for block in blocks:

41 for column in block.columns:

42 column.in_block = "+"

43 return "".join(column.in_block for column in alignment.columns)

45 def main():

46 alignment = dna.Alignment().append_file(infile)

47 assert len(alignment.sequences) == 2, "Input must have TWO sequences!"

48 alignment.realign(Left())

49 alignment.realign(Needle())

50 blocks = find_runs(alignment)

52 for n, block in enumerate(blocks, 1):

53 block.to_file(open("block_%02d.fasta" % n, "w"))

55 alignment.to_file(outfile)

56 FastaFile(outfile).write_string(

57 blocks_markup(alignment, blocks),

58 "markup",

59 "In run with window %s and threshold %s" % (width, threshold)

60 )

62 try:

63 main()

64 except Exception, e:

65 print "An error has occured:", e