DNA_ALPHABET does not include IUPAC symbols

[ArtPoon]

OpenRDP/openrdp/__init__.py

Line 29 in d5cb639

DNA_ALPHABET = ['A', 'T', 'G', 'C', '-', '*']

This can cause problems because sequences with mixtures (e.g., N or R) will cause the program to exit!

OpenRDP/openrdp/__init__.py

Lines 315 to 316 in d5cb639

	if not valid_alignment(aln) and not valid_chars(aln):
	sys.exit(1)

[ArtPoon]

Things might be more complicated than this - do any of the RDP methods reject sequences with mixtures?

[darrenpmartin]

If you're going tp do something to handle symbols other than ATGC- you should probably write separate versions of the main triplet scanning loops to do this - handling these in the inner loop will be very expensive.

…

On Wed, Feb 8, 2023 at 5:27 AM Art Poon ***@***.***> wrote: Things might be more complicated than this - do any of the RDP methods reject sequences with mixtures? — Reply to this email directly, view it on GitHub <#46 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ADEJ3TVOGK2B7FU4DKRMCZDWWMHDTANCNFSM6AAAAAAUUXQCPQ> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[ArtPoon]

Should be a simple modification of this function:

OpenRDP/openrdp/common.py

Lines 46 to 52 in d5cb639

	def percent_diff(s1, s2):
	s1_valid = (s1 == 'A') | (s1 == 'T') | (s1 == 'G') | (s1 == 'C')
	s2_valid = (s2 == 'A') | (s2 == 'T') | (s2 == 'G') | (s2 == 'C')
	valid = s1_valid & s2_valid
	diffs = np.sum(s1[valid] != s2[valid])
	num_valid = valid.sum()
	return diffs.sum() / num_valid if num_valid else 0

Have not yet found where characters - and * are being handled.

[ArtPoon]

Mixtures can be resolved or averaged at a pre-processing step, before looping over triplets.

[ArtPoon]

I just noticed that the test statement in __init__.py uses an and instead of an or, so it is letting through sequences that contain invalid characters anyhow!