in reply to Hamming Distance Between 2 Strings - Fast(est) Way?
#!/usr/bin/perl -slw
use strict;
my $s1 = 'AAAAA';
my $s2 = 'ATCAA';
my $s3 = 'AAAAA';
print "$s1:$s2 hd:", hd( $s1, $s2 ); # will give value 2
print "$s1:$s3 hd:", hd( $s1, $s3 ); # will give value 0
sub hd{ length( $_[ 0 ] ) - ( ( $_[ 0 ] ^ $_[ 1 ] ) =~ tr[\0][\0] ) }
__END__
P:\test>500235
AAAAA:ATCAA hd:2
AAAAA:AAAAA hd:0
Examine what is said, not who speaks -- Silence betokens consent -- Love the truth but pardon error.
Lingua non convalesco, consenesco et abolesco. -- Rule 1 has a caveat! -- Who broke the cabal?
"Science is about questioning the status quo. Questioning authority".
The "good enough" maybe good enough for the now, and perfection maybe unobtainable, but that should not preclude us from striving for perfection, when time, circumstance or desire allow.
Re^2: Hamming Distance Between 2 Strings - Fast(est) Way?
by Roy Johnson (Monsignor) on Oct 14, 2005 at 17:52 UTC
|
Or a little more plainly (update to explain: using the names the OP did, so he can see what corresponds):
sub hd {
my ($k, $l) = @_;
my $diff = $k ^ $l;
my $num_mismatch = $diff =~ tr/\0//c;
}
Caution: Contents may have been coded under pressure.
| [reply] [d/l] |
|
I'm not sure that $k, $l are any more meaningful than $_[0], $_[1] and you pay a small performance penalty for obtaining those names.
Does $diff really capture what it is naming? I've used $xor or $mask for that in the past, but I wonder if it isn't best left unnamed.
I don't see any purpose in naming the return value, just to return it. Better to say return <EXPR>; and give the function a proper name like hammingDistance() and allow it to name the return value.
While using tr/\0//c produces the same result as length - tr/\0/\0/, the major difference is that the former modifies the string being inspected, deleting the chars counted, where the latter does not. It just counts.
On the short strings in the OP this is insignificant. But the technique works for any length strings, and as DNA strings can get very large indeed, the difference then becomes very significant.
Overall, I prefer my version to yours, but each to their preference :)
Examine what is said, not who speaks -- Silence betokens consent -- Love the truth but pardon error.
Lingua non convalesco, consenesco et abolesco. -- Rule 1 has a caveat! -- Who broke the cabal?
"Science is about questioning the status quo. Questioning authority".
The "good enough" maybe good enough for the now, and perfection maybe unobtainable, but that should not preclude us from striving for perfection, when time, circumstance or desire allow.
| [reply] [d/l] [select] |
|
| [reply] [d/l] [select] |
|
Re^2: Hamming Distance Between 2 Strings - Fast(est) Way?
by b_hall (Scribe) on Oct 14, 2005 at 17:35 UTC
|
Depending on the length and number of your nucleic acids, and the number of comparisons you are making it might be better to store and compare bit strings generated with the vec command rather than character strings. You could convert each character in the original string with a for loop and a hash table ie something like
my $char_string = "AAATAT";
my $bit_string = "";
my %lookup_table = ("A" => "00", "T" => "01", "C" => "10", "G" => "11
+");
my @char_array=(split //, $char_string);
for (my $i=0; $i <= $#char_array ; $i++) {
vec($bit_string,$i,2) = $lookup_table{$char_array[$i]};
}
| [reply] [d/l] |
|
| [reply] |
|
I don't think that this would be a suitable solution for every problem of this type, but I think how useful it is will depend on the precise problem being addressed. This is not my field, but here are some thoughts on your reply.
One situation where you would want to compare multiple long strings would be if you have a set of sequenced functional mutants which you needed to compare to the wildtype, and to one another to see if you can cluster them and obtain insight into function. If there were enough mutants, this may promote the vec method. If you are looking at non standard bases (perhaps from damage or being from an unusual archaea) this solution might not be suitable, but my understanding was that usually non standard bases replace the common ones, not add to the list so you would still keep the total to four (one example that comes to mind is in RNA U, uracil, replaces T, thymine). Unusual bases would also only appear in a fraction of problems, so it might well not be a problem here.
What I don't know is how large the problem set would have to be to offset the encoding time. If somebody knows how to find that, that would help determine if the solution is appropriate for a given problem. Unfortunately I don't.
| [reply] |
|
|
|