the clearest and most efficient way to read FASTA file:
fh=open('seqA.fa')
name=fh.readline()[1:-1] #between '>' and '\n' is name of the sequence.
sequence=''
for line in fh:
sequence+=fh.readline()[:-1]
fh.close()
It is better to use readline() than read(). read() will read the whole file at once and store it as a string. This i a potential problem if there is not enough memory available to store the file's content. (Genome sequence files are usually in MB or GB size.)
Some modules in Biopython package allow us to do more complicated operations on sequence file if we need more information than the plain sequences, like the name, id, description, and cross references and annotations.
Here are some important modules:
SeqRecord
SeqIO , which is common interface to input and output sequence file formats.
fh=open('seqA.fa')
for record in SeqIO.parse(fh,'fasta'):
id=record.id
seq=record.seq
print id,'size:',len(seq)
fh.close()
Monday, November 14, 2011
Data types
New in Python 3: Dictionary
In python 3, keys(),values() and items() methods don't return a list as in previous versions. They return an iterable object. i.e. IUPAC={'A': 'Ala', 'C': 'Cys','E': 'Glu'}
>IUPAC.values() # under python 3.0, you will get:
<dict_values object at 0xb7x8s7s2>
>for x in IUPAC.values():
print x
Ala
Cys
Glu
This is more memory efficient than Python 2.x, since Python doesn't store the whole list in memory.
If you need to get the list anyway, use list()
>list(IUPAC.values())
['Ala', 'Cys', 'Glu']
In python 3, keys(),values() and items() methods don't return a list as in previous versions. They return an iterable object. i.e. IUPAC={'A': 'Ala', 'C': 'Cys','E': 'Glu'}
>IUPAC.values() # under python 3.0, you will get:
<dict_values object at 0xb7x8s7s2>
>for x in IUPAC.values():
print x
Ala
Cys
Glu
This is more memory efficient than Python 2.x, since Python doesn't store the whole list in memory.
If you need to get the list anyway, use list()
>list(IUPAC.values())
['Ala', 'Cys', 'Glu']
Subscribe to:
Posts (Atom)