Conservation scores at each location of few motifs in MSA
Entering edit mode
8 days ago

I am currently working on a bioinformatics problem where I need to lookup and count the location and count of occurences of 4000-ish 5 character long patterns in each sequence of a fasta file of 700GB. This fasta file has aligned sequences of one virus. To start with, I want to just look up 6 patterns.

The characters are limited to ATGC for the patterns and the fasta.

The input to the code is just the patterns I want to look up and the aligned sequences of the virus.

The output should be a file which is of a matrix format with #rows = num of patterns, # columns = #num of characters in any given fasta sequence (should be same for all sequences in the fasta as it is an aligned file). Each value of the matrix should be the % of sequences that have that pattern starting in that coordinate.

For example- input:

seq1 - ATGCCAT
seq2 - ATG-GAT
seq3 - GGAATTC
seq4 - -ATATCC
pattern #1 - AT
pattern #2 - GC


- = any char which has not been sequenced/ not in ATGC total number of sequences at any particular coordinate should be = number of seq - num of invalid sequences at that coordinate

output matrix:

Pattern 1 2 3 4 5 6
AT 2/3 1/4 0 2/3 0 2/4
GC 0 0 1/3 0 0 0

let me explain the value 2/3 in position 1*1 of the matrix - in the first coordinate of 2 out of the 4 sequences, the pattern AT starts. 1 sequence has an invalid - character and hence will be excluded from the total number of considered sequences.

I am a newbie to bioinformatics applications and hence am unaare if this is common problem but I am yet to discover a straghtforard solution through my research so far. I recently came across the aho-corasick algortihm for multi-pattern matching with strings through large text files as a very memory efficient and quick means of solving such issues. I am interested in scaling this search operation to many more patterns in the future as well. I am looking for advice on the following:

  • Are there any well- established tool for such and operation already? if so, what are they?
  • Am I right in attempting to use aho-corasick for this application?
  • How can I efficiently apply this algorithm to a largned fasta with aligned sequenced? can I read the fasta in chunks to make it quicker?
  • Should I code this in something other than Python? (I have been using Python so far but I am open to others as well).

Open to any suggestions as I am new to big data processing as well

conservation-scores FASTA MSA • 190 views
Entering edit mode

Login before adding your answer.

Traffic: 2109 users visited in the last hour
Help About
Access RSS

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6