I have a text file like this example:

>chr12:86512-86521

CGGCCAAAG

>chr16:96990-96999

CTTTCATTT

>chr16:97016-97025

TTTTGATTA

>chr16:97068-97077

ATTTAGGGA

This file is divided into different parts, and every part has 2 lines. The line which starts with > is ID and the 2nd line is a sequence of letters and the letters are A, T, C or G and also the length of each sequence is 9 so, for every sequence of letters there are 9 positions. I want to get the frequency of the 4 mentioned letters in every position (we have 9 positions).

Here is the expected output for the small example:

one = {'T': 1, 'A': 1, 'C': 2, 'G': 0}

two = {'T': 3, 'A': 0, 'C': 0, 'G': 1}

three = {'T': 3, 'A': 0, 'C': 0, 'G': 1}

four = {'T': 3, 'A': 0, 'C': 1, 'G': 0}

five = {'T': 0, 'A': 1, 'C': 2, 'G': 1}

six  = {'T': 0, 'A': 3, 'C': 0, 'G': 1}

seven = {'T': 2, 'A': 1, 'C': 0, 'G': 1}

eight = {'T': 2, 'A': 1, 'C': 0, 'G': 1}

nine = ({'T': 1, 'A': 2, 'C': 0, 'G': 1}

I am doing that in Python using the following command. This command has 3 steps. Steps 1 and 2 work fine, but would you help me to improve step 3, which made this pipeline so slow for big files?

Step 1: to parse the file into a comma-separated file

def fasta_to_textfile(filename, outfile): 

    with open(filename) as f, open(outfile, 'w') as outfile:

        header = sequence = None

        out = csv.writer(outfile, delimiter=',')

        for line in f:

            if line.startswith('>'):

                if header:

                    entry = header + [''.join(sequence)]

                    out.writerow(entry)

                header = line.strip('>n').split('|')

                sequence = 

            else:

                sequence.append(line.strip())

        if header:

            entry = header + [''.join(sequence)]

            out.writerow(entry)

Step 2: comma-separated file to a Python dictionary

def file_to_dict(filename):

    f = open(filename, 'r')

    answer = {}

    for line in f:

        k, v = line.strip().split(',')

        answer[k.strip()] = v.strip()

    return answer

To print functions from steps 1 and 2:

a = fasta_to_textfile('infile.txt', 'out.txt')

d = file_to_dict('out.txt')

Step 3: to get the frequency

one=

two=

three=

four=

five=

six=

seven=

eight=

nine=

mylist = d.values()

for seq in mylist:

    one.append(seq[0])

    two.append(seq[1])

    se.append(seq[2])

    four.append(seq[3])

    five.append(seq[4])

    six.append(seq[5])

    seven.append(seq[6])

    eight.append(seq[7])

    nine.append(seq[8])



from collections import Counter

one=Counter(one)    

two=Counter(two)

three=Counter(three)

four=Counter(four)

five=Counter(five)

You forgot to add import csv and from collections import Counter. Probably missed it while copy pasting. Also, your = signs are inconsistent in step 3. Try to follow PEP8. Also, a is useless in this line:

a = fasta_to_textfile('infile.txt', 'out.txt')

Since you've programmed a void function, a = None because it returns nothing.

Is the conversion to the CSV file really necessary? This would be an example of the pipeline:

1. Read the file.


2. Extract the sequence and load it into a N*9 table, where N is the number of sequences


3. Swap the rows and columns (numpy can help you out here)


4. A simple for loop that uses the Counter function on each row (but really column), refactored into less lines. Sadly I don't have time right to rewrite bits of your code right now.

One last thing - are you sure your example is correct? I tried loading it but got: ValueError: not enough values to unpack (expected 2, got 1)...