Working algos for Biological Data: Simple to Complex problems

''' Code description: Calculate Cys-richness of a protein with criteria set as: >=4 'C's over the length of protein AND >=5% total cysteine content Function: Take input sequence => Count number of 'C's & length => Calculate percentage Output True or False if criteria is met ''' from Bio import SeqIO def Cys_rich (record_seq): C_count = record_seq . count( 'C' ) seq_len = len (record_seq) Cys_perc = float (C_count) / float (seq_len) * 100 if C_count >= 4.0 and Cys_perc >= 5.0 : return 'Cys-rich' else : return 'No' CysRichSeq = [] for record in SeqIO . parse( 'filename.fasta' , 'fasta' ): if Cys_rich(record . seq) == 'Cys-rich' : CysRichSeq . append(record) SeqIO . write(CysRichSeq, 'Cys-rich_sequences.fasta' , 'fasta' ) print 'Cys-rich sequences written to file..'

A simple pandas solution to map multiple annotations for a protein.  A protein or gene file will have annotations curated by different methods. Most frequently, biologists will encounter more than one annotation for a single protein. It is a task in itself to pick the right annotation.  One of the simple ways is to consolidate them and pick the right ones after enough evidence is known.  Coding i n matlab or other languages may require more number of lines to achieve the same output. Here, a simple 'groupby' of pandas can do produce the same outputin seconds ! # Map multiple annotation for protein ID and join with a delimiter ''' Sample input A_xx Annotation1 A_xx Annotation2 B_xx Annotation1 Sample output A_xx Annotation1, Annotation2 B_xx Annotation1 ''' import pandas as pd data = pd . read_csv( 'input_file.txt' , delimiter = ' \t ' ) dfc = data . groupby([ 'PROTID' ])[ 'ANNOTATION...

#Simple code to find the occurrence of list of search terms in single line in a huge file. Search_terms = [ 'A' , 'B' , 'C' , 'D' ] with open ( 'BigFile.txt' , 'r' ) as infile : entries = infile . read () each_line = entries . splitlines () new_list = [] for ix , row in enumerate ( each_line ): element = row . split ( "\t" ) if set ( Search_terms ) == set ( element ): #Use <= if you want no strict option new_list . append ( element ) print ix + 1 #List the matching line number ! thefile = open ( 'Output.txt' , 'w' ) for item in new_list : print >> thefile , item

Few programs require certain versions of python. Especially, when you do not have root permission in your unix machine, you can still install a python version in your /home/usr directory. You can run it in parallel to in-built version. Follow these steps: 1. Unpack with  tar -xvf python-x.x.tar 2. ./configure 3. make altinstall prefix=~ exec-prefix=~ (./lib and ./bin will be created in your home directory(~)) 4. create alias for python-x.x => i.e cd ~/bin/ > ln -s python-x.x python ! If you create this outside bin directory, you will get a warning (ln: symbolic link already exists !! remember it refers to inbuilt python) 5. Complete alias creation by editing: vi ~/.bashrc with alias python='~/bin/python' Now any software you try to install with python setup.py it will access the version you installed in your home ! You can still install softwares using your version by skipping steps 4 & 5. But, everytime you need to keep specifying the export PATH ...

Simple, intutive & self-explanatory code to copy or remove columns in tab delimited text files. Remember AWK or Sed one liners are very handy too. But sometimes, if there are space or inconsistencies in file they may fail. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 __author__ = 'Arun Prasanna' ''' Remove first column of the txt file ''' with open ( 'Inputfile.txt' , 'r' ) as infile: entries = infile . read() each_line = entries . splitlines() new_list = [] for row in each_line: element = row . split( " \t " ) ele_size = len (element) for i in range ( 1 , ele_size): tmp = element[i] new_list . append(tmp) new_list . append( ' \t ' ) new_list . append( ' \n ' ) f = open ( 'Output.txt' , 'w' ) out = f . writelines(new_list) f . close() print...

Python code to count the number of elements (genes | proteins | genus ...) in each row in a non-homogenous cluster files. Example Input: g1  g2   g3  g4 g2 g4  g6  g7 Example Output: 4 1 3 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 __author__ = 'Arun Prasanna' ''' Python code to count number of elements in non-homogenous text file. Small, simple & self explanatory code ! ''' with open ( 'Input.txt' , 'r' ) as infile : entries = infile . read () each_line = entries . splitlines () new_list = [] for row in each_line : element = row . split ( "\t" ) ele_size = len ( element ) new_list . append ( str ( ele_size )) new_list . append ( '\n' ) f = open ( 'Count_EachRowElements.txt' , 'w' ) out = f . writelines ( new_list ) f . close () print "Program complete"

Code is useful to refine clustering data with IDs tagged with genus name. The output can be used to count protein copy numbers and etc., to create phyletic matrices or copy number matrices. Example Input: 123_g1   NID_g2   4567_g3   xx_g4    012_g6 NID_g10  ACC_g4 Example Output: g1    g2 g3    g4    g6 g10  g4 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 __author__ = 'Arun Prasanna' ''' Program to read the input file with 'genus_number' format names and convert that to => Each lines can have different number of elements (non-homogenous data) 1. 'genus_ID' to 'genus' format => [0] in split 2. 'ID_genus' to 'genus' format => [1] in split ''' with open ( 'Input_file.txt' , 'r' ) as infile : entries = infile . read () each_line = entries . splitlines () new_list ...

Convert the binary matrix to fasta format with this simple code in Python !. Recommended for larger file sizes ! Sample Input: Sp1 1 1 1 1 Sp2 0 0 0 0 Sp3 1 1 0 0 Sample Output: >Sp1 ['1', '1', '1', '1'] >Sp2 ['0', '0', '0', '0'] >Sp3 ['1', '1', '0', '0'] Workaround the output file in any text file to remove [,'  ] to generate final output as: >Sp1 1111 >Sp2 0000 >Sp3 1100 __author__ = 'Arun Prasanna' ''' This is a simple python code to convert binary matrix into fasta format. There  are many public softwares available. But each one has size limits (<=2MB).  This code processes 61 x 51000 character matrix in less than 10 seconds in  python 2.7! ''' with open ( 'infile_matrix-cp' , 'r' ) as infile : entries = infile . read () . strip () each_line = entries . ...