What is happening here


#1

i have some perl code

while ($high - $low > 0.001) { # precision
# calculate the sum of all normalized scores
my $sum = Pn * Pn * exp($lambda * $match) * 4
+ Pn * Pn * exp($lambda * $mismatch) * 12;
# refine guess at lambda
if ($sum > 1) {
$high = $lambda;
$lambda = ($lambda + $low)/2;
}
else {
$low = $lambda;
$lambda = ($lambda + $high)/2;
}
}

trying to rewrite this in ruby
while (high - low > 0.001)
#calculate the sum of all normalized scores
sum = Pn * Pn **(lambda * match) * 4 + Pn * Pn **(lambda * mismatch) *
12

#refine guess at lambda
if sum > 1
high = lambda
lambda = (lambda + low)/2
else
low = lambda
lambda = (lambda+high)/2
end
end

Is this the right way?

i seem to get different results when i execute the ruby version.
any ideas?


#2

George G. wrote:

sum = Pn * Pn **(lambda * match) * 4 + Pn * Pn **(lambda * mismatch) *
12

It seems like you’re using the ** operator wrong. You should use Ruby’s
Math.exp instead (which is the direct equivalent to Perl’s exp).

-Matthias


#3

It seems like you’re using the ** operator wrong. You should use Ruby’s
Math.exp instead (which is the direct equivalent to Perl’s exp).

-Matthias

Thank you Mat
i have changed that. here is the full script in ruby:
Pn = 0.25 #probabilty of any nuclotide

match = 10
mismatch = -10

expected_score = match * 0.25 + mismatch * 0.75
if match <=0 or expected_score >= 0
puts “illegal scores”

quit

end

#calculate lambda
#initial estimates
lambda = 1
high = 2
low = 0

while (high - low > 0.001) do
#calculate the sum of all normalized scores
sum = Pn * Pn * Math.exp(lambda * match) * 4 + Pn * Pn *
Math.exp(lambda * mismatch) * 12

#refine guess at lambda
if sum > 1
high = lambda
lambda = (lambda + low)/2
else
low = lambda
lambda = (lambda+high)/2
end
end

#compute target frequencies and H
target_id = Pn * Pn * Math.exp(lambda * match) * 4
h = lambda * match * target_id + lambda * mismatch * (1 - target_id)

puts “expected_score #{expected_score}\n”
puts “lambda: #{lambda} nats, #{lambda/Math.log(2)} bits”
puts “H: #{h} nats , #{h/Math.log(2)} bits”
puts “%ID: #{target_id * 100}”

and the original perl script is here:
#!/usr/bin/perl -w use strict;
use constant Pn => 0.25; # probability of any nucleotide

die “usage: $0 \n” unless @ARGV == 2;

my ($match, $mismatch) = @ARGV;

my $expected_score = $match * 0.25 + $mismatch * 0.75;
die “illegal scores\n” if $match <= 0 or $expected_score >= 0;

calculate lambda

my ($lambda, $high, $low) = (1, 2, 0); # initial estimates

while ($high - $low > 0.001) { # precision
# calculate the sum of all normalized scores
my $sum = Pn * Pn * exp($lambda * $match) * 4
+ Pn * Pn * exp($lambda * $mismatch) * 12;
# refine guess at lambda
if ($sum > 1) {
$high = $lambda;
$lambda = ($lambda + $low)/2;
}
else {
$low = $lambda;
$lambda = ($lambda + $high)/2;
}
}

compute target frequency and H

my $targetID = Pn * Pn * exp($lambda * $match) * 4;
my $H = $lambda * $match * $targetID
+ $lambda * $mismatch * (1 -$targetID);

output

print “expscore: $expected_score\n”;
print “lambda: $lambda nats (”, $lambda/log(2), " bits)\n";
print “H: $H nats (”, $H/log(2), " bits)\n";
print "%ID: ", $targetID * 100, “\n”;

The ruby version seems to get into an infinite loop. while the perl
version executes nicely and with the expected results.


#4

Quoting George G. removed_email_address@domain.invalid:

}
else {
    $low = $lambda;
    $lambda = ($lambda + $high)/2;
}

}

Pn * exp($lambda * $match) is Pn * (Math::E ** (lambda + match)), not
Pn ** (lambda + match).

WBR, Peter Z.


#5

define lambda, high and low as floats, lambda = 1.0, high=2.0, low = 0.0

Am Mittwoch, den 08.04.2009, 07:57 -0500 schrieb George G.:


#6

Fritz Heinrichmeyer wrote:

define lambda, high and low as floats, lambda = 1.0, high=2.0, low = 0.0

Thank you so much everyone for your responses, it works nicely now!

GG


#7

in ruby lambda, high and low are integers and divide by two is also
integer arithmetic. Perl uses float for scalars by default.

Am Mittwoch, den 08.04.2009, 07:57 -0500 schrieb George G.: