Context:
This script reads in a fasta file, counts the number '-' characters, alphabetic characters and total characters for each sequence (sequence descriptions are on lines marked with '>'), and then prints out a processed file with columns where 90% (or more) of '-' characters are removed.
I've tested it with smaller alignments and it worked fine, but it's currently been running for around 4 hours on one alignment (220Mb), which seems slow and probably inefficient.
The script:
import sys
from Bio import AlignIO
alignment = AlignIO.read(open(sys.argv[1]), "fasta")
length = alignment.get_alignment_length()
tofilter = []
for i in range(length):
counterofgaps=0
counterofsites=0
counteroftotal=0
for record in alignment:
for site in record.seq[i]:
if site == '-':
counterofgaps +=1
counteroftotal +=1
elif site.isalpha():
counterofsites +=1
counteroftotal +=1
Fractionofgaps = counterofgaps/counteroftotal * 100
if Fractionofgaps >=90:
tofilter.append(i)
for record in alignment:
print ('>' + record.id)
for i in range(length):
if i not in tofilter:
print (record.seq[i], end ='')
print()
SeqIO.parse
instead help? \$\endgroup\$-
char. Potentially, all records that fall into crucial condition ">= 90
% of-
chars " can have those-
chars (gaps) as the rightmost sequence and alpha chars - as leftmost, at the very start. But in such case - all alpha letters would be skipped. Why is that logic correct? \$\endgroup\$