Based on your work, you can process fasta file in anyway. But through this post, I am going to present a simple example to find the average sequence length from a fasta file . Reading a fasta file is already discussed.
Program to find mean value from given fasta file
In order to find mean value, we don’t need header information. Therefore a simple python list is enough to store base pair / amino acid count of each sequence. Total number of items in the list represent total number sequence present in the file.
If we declare a list ‘fasta’ to store base pair / amino acid count from each sequence, then length of ‘fasta’ return the number of sequence present in the file.
fasta = [23,36,54] # sample list after loading file
totalSequence = len(fasta)
print(totalBase) # output : 3
Since ‘fasta’ list store base pair / amino acid count from each sequence, passing ‘fasta’ into python3 sum() function return total base pair / amino acids present in the file.
Fasta file is a text-based file used to store nucleotide or protein sequences.Inside the file, each sequence begins with a single line description. The description line starts with a greater than (‘>‘) symbol. Description is followed by lines containing sequence data.
Example:
Reading a file using python
Python programming language is simple to use. Using python3, we can simply read a text file using following few lines of code:
file = open("location/to/file.txt",'r')
for line in file:
print(line)
file.close()
The code is simple. Read file inside the variable ‘file’, Start loop for each ‘line’ inside the ‘file’, print that ‘line’, Close the file.
This can also be done in following way:
with open("location/to/file.txt",'r') as file:
for line in file:
print(line)
The advantage of above code is that, it automatically close the file once the processing is done.
Reading fasta file using python
In order to do this, we need a python dictionary to store our sequence data. A python dictionary is an unordered collection of items, where each item is accessed by referring to its key name.
For fasta file, we are going to specify description of sequence as key and sequence as item.
We know that every description start with a ‘>’ symbol. Thus, by checking if a line starts with ‘>’ symbol, we can easily identify description line and sequence line.
When a line is found to be start with ‘>’ symbol, the whole line is considered as a key inside the dictionary. Until a new description line is found, all those sequence line followed by current description line is appended to its value.
While processing, the hidden next line character ‘\n’ in the end of each line is removed. This is accomplished simply by cropping the end character of each line.
By combining all these logic, we can read a fasta file into a dictionary using following code:
inputfilename = 'test.fasta'
output_dictionary = {}
with open(inputfilename, 'r') as file:
current_header = ''
for line in file:
if line.startswith('>'):
current_header = line[1:][:-1]
output_dictionary[current_header] = ""
else:
if not(current_header == ''):
if line.endswith('\n'):line = line[:-1]
output_dictionary[current_header] += line
In order to avoid the loss of last sequence when the loop ends, an extra line of code is required outside the loop, output_dictionary[current_header] += line make sure we never miss our last sequence.
