Unable to parse just sequences from FASTA file

Question

How can I remove ids like '>gi|2765658|emb|Z78533.1|CIZ78533 C.irapeanum 5.8S rRNA gene and ITS1 and ITS2 DNA\n' from sequences?

I have this code:

with open('sequence.fasta', 'r') as f :
    while True:
        line1=f.readline()
        line2=f.readline()
        line3=f.readline()
        if not line3:
            break
        fct([line1[i:i+100] for i in range(0, len(line1), 100)])
        fct([line2[i:i+100] for i in range(0, len(line2), 100)])
        fct([line3[i:i+100] for i in range(0, len(line3), 100)])

Output:

['>gi|2765658|emb|Z78533.1|CIZ78533 C.irapeanum 5.8S rRNA gene and ITS1 and ITS2 DNA\n']
['CGTAACAAGGTTTCCGTAGGTGAACCTGCGGAAGGATCATTGATGAGACCGTGGAATAAACGATCGAGTG\n']
['AATCCGGAGGACCGGTGTACTCAGCTCACCGGGGGCATTGCTCCCGTGGTGACCCTGATTTGTTGTTGGG\n']
['CCGCCTCGGGAGCGTCCATGGCGGGTTTGAACCTCTAGCCCGGCGCAGTTTGGGCGCCAAGCCATATGAA\n']
['AGCATCACCGGCGAATGGCATTGTCTTCCCCAAAACCCGGAGCGGCGGCGTGCTGTCGCGTGCCCAATGA\n']
['ATTTTGATGACTCTCGCAAACGGGAATCTTGGCTCTTTGCATCGGATGGAAGGACGCAGCGAAATGCGAT\n']
['AAGTGGTGTGAATTGCAAGATCCCGTGAACCATCGAGTCTTTTGAACGCAAGTTGCGCCCGAGGCCATCA\n']
['GGCTAAGGGCACGCCTGCTTGGGCGTCGCGCTTCGTCTCTCTCCTGCCAATGCTTGCCCGGCATACAGCC\n']
['AGGCCGGCGTGGTGCGGATGTGAAAGATTGGCCCCTTGTGCCTAGGTGCGGCGGGTCCAAGAGCTGGTGT\n']
['TTTGATGGCCCGGAACCCGGCAAGAGGTGGACGGATGCTGGCAGCAGCTGCCGTGCGAATCCCCCATGTT\n']
['GTCGTGCTTGTCGGACAGGCAGGAGAACCCTTCCGAACCCCAATGGAGGGCGGTTGACCGCCATTCGGAT\n']
['GTGACCCCAGGTCAGGCGGGGGCACCCGCTGAGTTTACGC\n']
['\n']
...

My function is:

def fct(input_string):
    code={"a":0,"c":1,"g":2,"t":3}
    p=[code[i] for i in input_string]
    n=len(input_string)
    c=0

    for i, n in enumerate(range(n, 0, -1)):
        c +=p[i]*(4**(n-1))
        return c+1

fct() returns an integer from a string. For example, ACT gives 8 i.e.: my function must take as input string sequences contain just the following bases A,C,G,T

But when I use my function it gives:

KeyError: '>gi|2765658|emb|Z78533.1|CIZ78533 C.irapeanum 5.8S rRNA gene and ITS1 and ITS2 DNA\n'

I try to remove ids by stripping lines start with > and writing the rest in text file so, my text file output.txt contains just sequences without ids, but when I use my function fct I found the same error:

KeyError: 'CGTAACAAGGTTTCCGTAGGTGAACCTGCGGAAGGATCATTGATGAGACCGTGGAATAAACGATCGAGTG\n'

What can I do?

On an unrelated note, you should consider giving [PEP8](http://www.python.org/dev/peps/pep-0008/) a read to clean up your Python formatting. Most importantly, use spaces, not tabs! — David Cain, Jul 30 '13 at 00:11
Please don't delete your comments and questions - these sorts of questions are a helpful resource to others with similar problems. — David Cain, Aug 16 '13 at 19:48

David Cain · Answer 1 · 2013-07-30T00:27:37.863

I see two major problems in your code: You're having problems parsing FASTA sequences, and your function is not properly iterating over each sequence.

Parsing FASTA data

Might I suggest using the excellent Biopython package? It has excellent FASTA support (reading and writing) built in (see Sequences in the Tutorial).

To parse sequences from a FASTA file:

for seq_record in SeqIO.parse("seqs.fasta", "fasta"):
    print record.description  # gi|2765658|emb|Z78533.1...
    print record.seq  # a Seq object, call str() to get a simple string

>>> print record.id
'gi|2765658|emb|Z78533.1|CIZ78533'

>>> print record.description
'gi|2765658|emb|Z78533.1|CIZ78533 C.irapeanum 5.8S rRNA gene and ITS1 and ITS2 DNA'

>>> print record.seq
Seq('CGTAACAAGGTTTCCGTAGGTGAACCTGCGGAAGGATCATTGATGAGACCGTGG...CGC', SingleLetterAlphabet())

>>> print str(record.seq)
'CGTAACAAGGTTTCCGTAGGTGAACCTGCGGAAGGATCATTGATGAGACC'  #(truncated)

Iterating over sequence data

In your code, you have a list of strings being passed to fct() (input_string is not actually a string, but a list of strings). The solution is just to build one input string, and iterate over that.

Other errors in `fct`:

You need to capitalize the keys to your dictionary: case matters
You should have the return statement after the for loop. Keeping it nested means c is returned immediately.
Why bother constructing p when you can just index into code when iterating over the sequence?
You write over the sequence's length (n) by using it in your for loop as a variable name

Modified code (with proper PEP 8 formatting), and variables renamed to be clearer what they mean (still have no idea what c is supposed to be):

from Bio import SeqIO


def dna_seq_score(dna_seq):
    nucleotide_code = {"A": 0, "C": 1, "G": 2, "T": 3}

    c = 0 
    for i, k in enumerate(range(len(dna_seq), 0, -1)):
        nucleotide = dna_seq[i]
        code_num = nucleotide_code[nucleotide]
        c += code_num * (4 ** (k - 1)) 
    return c + 1 


for record in SeqIO.parse("test.fasta", "fasta"):
    dna_seq_score(record.seq)

Actually, if you try running this code, you'll see a different error message. I've updated my answer to explain what's going on here. — David Cain, Jul 30 '13 at 00:19
Perhaps this question will be of use to you: http://stackoverflow.com/questions/312443/how-do-you-split-a-list-into-evenly-sized-chunks-in-python — David Cain, Jul 30 '13 at 00:41
Okay. Have you tried converting your code to work with a `Seq` object, or the sequence as a string? If you still have issues, post another question, and I'm sure you'll get help resolving them. — David Cain, Jul 30 '13 at 01:34
You'll need to be more specific. Either [edit](http://stackoverflow.com/posts/17936166/edit) your question with the exact code you're using, or [ask a new question](http://stackoverflow.com/questions/ask). — David Cain, Jul 30 '13 at 12:33

Unable to parse just sequences from FASTA file

1 Answers1

Parsing FASTA data

Iterating over sequence data

Other errors in `fct`:

Linked

Unable to parse just sequences from FASTA file

1 Answers1

Parsing FASTA data

Iterating over sequence data

Other errors in fct:

Linked

Other errors in `fct`: