Answers.Body: <p>I wonder if FIMO might be something that you could use.
  The online version doesn't have the capacity 70 deal with 700GB of data, but you could try it out with a subsample to see if the output is something that would fit your needs. <a href="https://meme-suite.org/meme/tools/fimo" rel="nofollow noreferrer">https://meme-suite.org/meme/tools/fimo</a>
  The MEME suite (including FIMO) is also available via conda, or Docker etc.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2024-05-30T12:36:05.987Z
Answers.LastActivityDate: 2024-05-30T12:36:05.987Z
Answers.UserDisplayName: Laura
Answers.Body: <p>Maybe a job for <code>seqkit locate</code>. Here <code>-m 0</code>means no mismatches allowed. I do not understand really your output format, so maybe start with seqkit and see whether you can tune it for your problem. It allows mismatches by some parameters.</p>
  <pre><code>$ cat foo.fa 
  &gt;entry1
  TACTACGATCACTGACGTAGCACATTGGATTTAGAGTACGATCGAT
  &gt;entry2
  TAGCACTATCGACGTAGCTACGTACGTACTGATCGATCTGACACGA
  &gt;entry3
  TGGTGTGTTGGTGTGTTGAAAAAAAAAAAAAAAAAAAAAAAAAAAA
  
  $ cat pattern.fa 
  &gt;pattern1
  TAGCAC
  
  $ ./seqkit locate -m 0 -f pattern.fa foo.fa 
  seqID   patternName pattern strand  start   end matched
  entry1  pattern1    TAGCAC  +   18  23  TAGCAC
  entry2  pattern1    TAGCAC  +   1   6   TAGCAC
  </code></pre>
  <p><a href="https://bioinf.shenwei.me/seqkit/usage/#locate" rel="nofollow noreferrer">https://bioinf.shenwei.me/seqkit/usage/#locate</a></p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2024-05-30T20:41:41.383Z
Answers.LastActivityDate: 2024-05-30T20:41:41.383Z
Answers.UserDisplayName: ATpoint
Body: <p>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.</p>
  <p>The <strong>characters are limited to ATGC</strong> for the patterns and the fasta.</p>
  <p>The <strong>input</strong> to the code is just the patterns I want to look up and the aligned sequences of the virus.</p>
  <p>The <strong>output</strong> 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.</p>
  <p><strong>For</strong> <strong>example</strong>-<br />
  <strong>input</strong>:</p>
  <p>seq1 - ATGCCAT<br />
  seq2 - ATG-GAT<br />
  seq3 - GGAATTC<br />
  seq4 - -ATATCC</p>
  <p>pattern #1 - AT<br />
  pattern #2 - GC</p>
  <p>NUMBER OF SEQUENCES = 4<br />
  '-' = 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</p>
  <p><strong>output</strong> matrix:</p>
  <p><strong>Pattern</strong> <strong>1 2 3 4 5 6</strong><br />
  <strong>AT</strong> 2/3 1/4 0 2/3 0 2/4<br />
  <strong>GC</strong> 0 0 1/3 0 0 0</p>
  <p><strong>let me explain the value 2/3</strong> in <strong>position 1*1 of the matrix</strong> - 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.</p>
  <p>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:</p>
  <ol>
  <li>Are there any well- established tool for such and operation already? if so, what are they?</li>
  <li>Am I right in attempting to use aho-corasick for this application?</li>
  <li>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?</li>
  <li>Should I code this in something other than Python? (I have been using Python so far but I am open to others as well).</li>
  </ol>
  <p>Open to any suggestions as I am new to big data processing as well</p>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-05-30T20:32:48.417Z
Comments.Text: This question needs more context about the overall project, because it seems like an odd task to carry out. What is the purpose of doing this? What will you do with these match results once the answer has been found?
Comments.UserDisplayName: gringer
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-05-30T21:04:34.110Z
Comments.Text: Wait ... **700GB** ... I have not read the question in detail, but I have produced solutions involving industrial data engineering (I work in industry), but you'd need to be an expert and it's an entire project. Thus it is a full scale parallelisation **OR** a single CPU with massive RAM (AWS stuff). I'll think if there's an easy way in ... BTW is 700GB zipped or unzipped (most fasta/fastq are zipped))?
Comments.UserDisplayName: M__
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-07T03:39:17.220Z
Comments.Text: The unzipped ones is 700GB
Comments.UserDisplayName: Swathi Subramanyan
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-07T03:39:50.323Z
Comments.Text: The purpose is to calculate conservation scores
Comments.UserDisplayName: Swathi Subramanyan
Tags: <python><fasta><algorithms><homework><multiple-sequence-alignment>
Title: Multi-pattern search in aligned sequences
ThreadText: # Multi-pattern search in aligned sequences
  
  (Asked by: username_0 on 2024-05-30)
  
  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
  
  NUMBER OF SEQUENCES = 4  
  '-' = 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:
  
  1. Are there any well- established tool for such and operation already? if so, what are they?
  2. Am I right in attempting to use aho-corasick for this application?
  3. 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?
  4. 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
  
  **Tags:** python, fasta, algorithms, homework, multiple-sequence-alignment
  
  > **Comments:**
  > **username_1**: This question needs more context about the overall project, because it seems like an odd task to carry out. What is the purpose of doing this? What will you do with these match results once the answer has been found?  
  > **username_2**: Wait ... **700GB** ... I have not read the question in detail, but I have produced solutions involving industrial data engineering (I work in industry), but you'd need to be an expert and it's an entire project. Thus it is a full scale parallelisation **OR** a single CPU with massive RAM (AWS stuff). I'll think if there's an easy way in ... BTW is 700GB zipped or unzipped (most fasta/fastq are zipped))?  
  > **username_0**: The unzipped ones is 700GB  
  > **username_0**: The purpose is to calculate conservation scores  
  
  ---
  
  ## 2 answers
  
  ### Answer by username_3 on 2024-05-30. Score: 1
  
  I wonder if FIMO might be something that you could use.
  The online version doesn't have the capacity 70 deal with 700GB of data, but you could try it out with a subsample to see if the output is something that would fit your needs. <https://meme-suite.org/meme/tools/fimo>
  The MEME suite (including FIMO) is also available via conda, or Docker etc.
  
  
  ---
  
  ### Answer by username_4 on 2024-05-30. Score: 1
  
  Maybe a job for `seqkit locate`. Here `-m 0`means no mismatches allowed. I do not understand really your output format, so maybe start with seqkit and see whether you can tune it for your problem. It allows mismatches by some parameters.
  
  ```
  $ cat foo.fa 
  >entry1
  TACTACGATCACTGACGTAGCACATTGGATTTAGAGTACGATCGAT
  >entry2
  TAGCACTATCGACGTAGCTACGTACGTACTGATCGATCTGACACGA
  >entry3
  TGGTGTGTTGGTGTGTTGAAAAAAAAAAAAAAAAAAAAAAAAAAAA
  
  $ cat pattern.fa 
  >pattern1
  TAGCAC
  
  $ ./seqkit locate -m 0 -f pattern.fa foo.fa 
  seqID   patternName pattern strand  start   end matched
  entry1  pattern1    TAGCAC  +   18  23  TAGCAC
  entry2  pattern1    TAGCAC  +   1   6   TAGCAC
  
  ```
  
  <https://bioinf.shenwei.me/seqkit/usage/#locate>
  
  
  ---
  
  