I have this function:
import re def check_and_clean_sequence(sequence, alphabet): """ Function to check and clean up all ambiguous bases in a sequence. Ambigous bases are bases that are not in the sequence alphabet, ie. 'ACGT' for DNA sequences. Inputs: sequence - a string representing a DNA sequence. alphabet - set/string representing the allowed characters in a sequence. Outputs: cleaned_sequence - cleaned sequence or a string. """ if set(sequence).issubset(alphabet): return sequence else: return cleaning_ambiguous_bases(sequence) def cleaning_ambiguous_bases(sequence): """ Function to clean up all ambiguous bases in a sequence. Ambiguous bases are bases that are not in the sequence alphabet, ie. 'ACGT' for DNA sequences. Inputs: sequence - a string representing a DNA sequence. Outputs: integer - a new clean up string representing a DNA sequence without any ambiguous characteres. """ # compile the regex with all ambiguous bases pat = re.compile(r'[NRYWXSKM]') # look for the ambiguous bases and replace by # nothing new_seq = re.sub(pat, '', sequence) return new_seq
The performance in jupyter notebook with timeit is around:
%timeit check_and_clean_sequence(seq, iupac_dna) 200 ms ± 436 µs per loop (mean ± std. dev. of 7 runs, 1 loop each)
This was calculated with a DNA sequence with 5539776 pb or 5539.776 kpb. The sequence can be downloaded from: E. coli genome
The function receive a big string and check if there are any characters that are not in the allowed alphabet (in the case the IUPAC DNA = 'ACGT'). If the sequence has they are cleaned up and the function returns a new sequence, otherwise return the original one.
Are there are tips or tricks to improve this time?
>>> import re >>> alphabet = 'ACGT' >>> toy_sequence = 'AGCTAGGGACATTGAGGACCACCACAXYXMAAAAAA' >>> check_and_clean_sequence(toy_sequence, alphabet) 'AGCTAGGGACATTGAGGACCACCACAAAAAAA'