How to register the ID and the count (or number) of occurrences in a sequence?

Viewed 194

I am very new to python so bear with me, please. I have been asked to create a program. I need to look for occurrences of a pattern (which has to be provided by the user using the keyboard) string inside each of the sequences and if there are occurrences, register the id and the number (count) of occurrences in the sequence.

The data looks like this:

id  sequence
1   MVLSEGEWAAVLHVWAKVEADVAAGHGQDILIRLFKS
2   MNIFEMLRIAAGLRLKIYKDTEAAGYYTIGIGHLLTKSPSL
3   MVLSEGEWQLVLHVWAKVEADVAGHGQDILIRLFKSH
4   MNIFEMLRAAEGAALRLKIYKAADTEGYYTIGIGHLLTKS
5   MVLSAAEGEWQLVLHVWAKVEADVAGHGQDILIRLFK

where the ids are the numbers and the sequence are the sequences below each number. The file is a matrix of (100437 x 2)

This is the code that I have so far:

import re

def proteins_pattern_count(pattern):
    with open("proteins.csv", 'r') as proteins:
        proteins = proteins.read()
        items = re.findall(pattern, proteins)
    return len(items)

# Reading the pattern to look for and forcing the input to change the pattern to capital letters.
pattern = input("Please type in the pattern you would like to look for: ").upper()

count = proteins_pattern_count(pattern)

print('The pattern {} appears {} times within the proteins file'.format(pattern, count))

Output that I get:

Please type in the pattern you would like to look for: AA
The pattern AA appears 173372 times within the proteins file

But what I really want is: for instance, if the pattern I am looking for is "AA" then I would like to see a table only with the ids of the sequences that actually have such pattern inside and the number (count) of occurrences in the sequence, like this:

id  count
1   2
2   2
4   3
5   1

I think this is very easy to do, but I'm completely new to Python.

Thanks for your support!!

2 Answers

First, I would highly recommend using a different variable name for your object returned by .readlines(), as currently you are overwriting your file pointer. While it isn't necessary for this code, it's good practice in general.

Second, to make things easier, you may wish to use csv.reader to split up your csv in a nice, easy fashion.

Here are three code snippets you can use:

import re, csv

def proteins_slow1(pattern):
    with open("proteins.csv", 'r') as fp:
        proteins = csv.reader(fp, delimiter=',')
        ids, counts = [], []
        for i,seq in proteins:
            count = len(re.findall(pattern, seq))
            if count != 0:
                ids.append(i)
                counts.append(count)
    return ids, counts

def proteins_slow2(pattern):
    with open("proteins.csv", 'r') as fp:
        proteins = csv.reader(fp, delimiter=',')
        struct = {}
        for i,seq in proteins:
            count = len(re.findall(pattern, seq))
            if count != 0:
                struct[i] = count
    return struct

def proteins_fast(pattern):
    with open("proteins.csv", 'r') as fp:
        proteins = csv.reader(fp, delimiter=',')
        struct = {}
        for i,seq in proteins:
            count = seq.count(pattern)
            if count != 0:
                struct[i] = count
    return struct

proteins_slow1 makes two lists, and appends them if the count is non-zero. The function ends by returning a tuple with the ids list and the counts list.

This is pretty much as fast as proteins_slow2, which makes a dictionary, and adds new entries as key-value pairs (id as the key, count as the value).

The quickest is to actually not use re but instead use the .count() method on the sequence string. This shaves off about 30-40% off the running time (which will become important if you're repeatedly looking at 100000+ rows).

>>> %timeit proteins_slow1('AA')
199 ms ± 12.1 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)
>>> %timeit proteins_slow2('AA')
187 ms ± 767 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)
>>> %timeit proteins_fast('AA')
119 ms ± 539 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)

(Timing tests were done by making a csv file from the following function:)

def writer():
    with open("proteins.csv", 'w', newline='\n', encoding='utf-8') as fp:
        proteins = csv.writer(fp, delimiter=',')

        for i in range(20000):
            proteins.writerow([str(0+i*5),'MVLSEGEWAAVLHVWAKVEADVAAGHGQDILIRLFKS'])
            proteins.writerow([str(1+i*5),'MNIFEMLRIAAGLRLKIYKDTEAAGYYTIGIGHLLTKSPSL'])
            proteins.writerow([str(2+i*5),'MVLSEGEWQLVLHVWAKVEADVAGHGQDILIRLFKSH'])
            proteins.writerow([str(3+i*5),'MNIFEMLRAAEGAALRLKIYKAADTEGYYTIGIGHLLTKS'])
            proteins.writerow([str(4+i*5),'MVLSAAEGEWQLVLHVWAKVEADVAGHGQDILIRLFK'])

Enjoy!

Something like this should work

res = {}
for protein in proteins.readlines():
    if pattern in protein:
        res[protein[0]] = protein.count(pattern)

print(res)
Related