There's more than one way to do things PerlMonks

### Re: Help with my Code

by poj (Prior)
 on Aug 10, 2013 at 13:23 UTC ( #1048901=note: print w/replies, xml ) Need Help??

in reply to Reaped: Help with my Code

Why on the first interation do you expect k = 44.530219 ? I calculate it to be 44.51357658.
poj

Replies are listed 'Best First'.
Re^2: Help with my Code
by BrowserUk (Pope) on Aug 10, 2013 at 14:04 UTC

Ditto. I suspect that the values in @frequency are wrong.

```#! perl -slw
use strict;
use List::Util qw[ sum ];

my @matrix =(
[0.54, 0.03, 0.03, 0.03, 0.08, 0.04, 0.05, 0.05, 0.02, 0.02, 0.03, 0.0
+4, 0.03, 0.01, 0.06, 0.11, 0.06, 0.01, 0.01,0.06],
[0.02, 0.57, 0.03, 0.01, 0.02, 0.06, 0.03, 0.01, 0.04, 0.01, 0.01, 0.1
+1, 0.02, 0.01, 0.01, 0.02, 0.02, 0.01, 0.01, 0.01],
[0.01, 0.02, 0.48, 0.06, 0.01, 0.03, 0.03, 0.03, 0.05, 0.01, 0.01, 0.0
+3, 0.01, 0.00, 0.01, 0.04, 0.03, 0.00, 0.02, 0.01],
[0.02, 0.01, 0.09, 0.57, 0.01, 0.04, 0.11, 0.02, 0.03, 0.00, 0.00, 0.0
+3, 0.01, 0.00, 0.02, 0.04, 0.02, 0.00, 0.01, 0.00],
[0.01, 0.00, 0.00, 0.00, 0.57, 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.0
+0, 0.00, 0.01, 0.00, 0.01, 0.01, 0.01, 0.01, 0.01],
[0.02, 0.04, 0.03, 0.02, 0.01, 0.42, 0.05, 0.01, 0.04, 0.01, 0.01, 0.0
+5, 0.02, 0.00, 0.01, 0.02, 0.02, 0.01, 0.01, 0.01],
[0.04, 0.04, 0.05, 0.13, 0.01, 0.11, 0.52, 0.02, 0.03, 0.01, 0.01, 0.0
+7, 0.01, 0.00, 0.02, 0.04, 0.04, 0.01, 0.01, 0.01],
[0.05, 0.02, 0.05, 0.03, 0.02, 0.02, 0.02, 0.77, 0.02, 0.00, 0.00, 0.0
+2, 0.01, 0.00, 0.01, 0.04, 0.01, 0.01, 0.00, 0.01],
[0.01, 0.01, 0.02, 0.01, 0.01, 0.02, 0.01, 0.00, 0.54, 0.00, 0.00, 0.0
+1, 0.01, 0.02, 0.01, 0.01, 0.01, 0.02, 0.04, 0.00],
[0.02, 0.01, 0.01, 0.00, 0.03, 0.01, 0.01, 0.00, 0.01, 0.50, 0.10, 0.0
+1, 0.09, 0.04, 0.01, 0.01, 0.03, 0.02, 0.01, 0.18],
[0.03, 0.02, 0.01, 0.01, 0.03, 0.03, 0.01, 0.00, 0.02, 0.13, 0.60, 0.0
+2, 0.19, 0.09, 0.02, 0.02, 0.02, 0.02, 0.03, 0.08],
[0.03, 0.14, 0.06, 0.03, 0.01, 0.09, 0.06, 0.02, 0.05, 0.01, 0.01, 0.5
+0, 0.02, 0.01, 0.02, 0.04, 0.03, 0.01, 0.01, 0.01],
[0.01, 0.01, 0.00, 0.00, 0.01, 0.02, 0.00, 0.00, 0.01, 0.03, 0.05, 0.0
+1, 0.42, 0.02, 0.00, 0.01, 0.01, 0.01, 0.01, 0.02],
[0.01, 0.01, 0.00, 0.00, 0.02, 0.01, 0.00, 0.00, 0.03, 0.02, 0.04, 0.0
+0, 0.04, 0.61, 0.00, 0.01, 0.01, 0.07, 0.14, 0.01],
[0.03, 0.01, 0.01, 0.01, 0.00, 0.01, 0.01, 0.01, 0.02, 0.00, 0.01, 0.0
+1, 0.01, 0.00, 0.72, 0.03, 0.01, 0.00, 0.00, 0.01],
[0.07, 0.03, 0.06, 0.04, 0.05, 0.04, 0.03, 0.03, 0.02, 0.01, 0.01, 0.0
+4, 0.02, 0.01, 0.04, 0.43, 0.09, 0.02, 0.01, 0.01],
[0.04, 0.02, 0.04, 0.02, 0.04, 0.03, 0.03, 0.01, 0.02, 0.02, 0.02, 0.0
+3, 0.03, 0.01, 0.02, 0.09, 0.52, 0.01, 0.01, 0.04],
[0.00, 0.00, 0.00, 0.00, 0.01, 0.00, 0.00, 0.00, 0.01, 0.00, 0.00, 0.0
+0, 0.00, 0.02, 0.00, 0.00, 0.00, 0.71, 0.02, 0.00],
[0.00, 0.01, 0.01, 0.00, 0.01, 0.01, 0.00, 0.00, 0.05, 0.01, 0.01, 0.0
+1, 0.01, 0.10, 0.00, 0.01, 0.01, 0.06, 0.61, 0.01],
[0.05, 0.01, 0.01, 0.01, 0.07, 0.02, 0.02, 0.01, 0.01, 0.20, 0.07, 0.0
+1, 0.06, 0.03, 0.01, 0.02, 0.06, 0.02, 0.02, 0.51]
);

my @freq = (
0.08477395, 0.05103334, 0.03837665, 0.05740129, 0.01256165,
0.03471746, 0.06883297, 0.07907659, 0.02077239, 0.06813598,
0.08906677, 0.06501537, 0.02318017, 0.03823936, 0.04036729,
0.05668318, 0.05721648, 0.00871765, 0.02785317, 0.07797831,
);

my @lnPAM = map[ map log( \$_ || 1e-308 ), @\$_ ], @matrix;
my \$k = 1;
my \$mutRate = sum map{ \$freq[ \$_ ] * ( 1 - \$matrix[ \$_ ][ \$_ ] ) } 0..
+19;
my \$iter = 0;

# The following loop goes until the difference reach the computer's
# precision.
while( abs( \$mutRate - 0.01) > 1e-16 ) {
++\$iter;
\$k *= \$mutRate * 100;
map \$_ /= ( \$mutRate * 100), @\$_ for @lnPAM;
@matrix = map[ map exp( \$_ ), @\$_ ], @lnPAM;
\$mutRate = sum map{ \$freq[ \$_ ] * ( 1 - \$matrix[ \$_ ][ \$_ ] ) } 0.
+.19;
printf "After %d iteration(s), k=%f, RateMutation=%f\n", \$iter, \$k
+, \$mutRate;
}

__END__
C:\test>junk
After 1 iteration(s), k=44.513577, RateMutation=0.013418
After 2 iteration(s), k=59.727332, RateMutation=0.010018
After 3 iteration(s), k=59.836899, RateMutation=0.010000
After 4 iteration(s), k=59.837487, RateMutation=0.010000
After 5 iteration(s), k=59.837490, RateMutation=0.010000
After 6 iteration(s), k=59.837490, RateMutation=0.010000
After 7 iteration(s), k=59.837490, RateMutation=0.010000

With the rise and rise of 'Social' network sites: 'Computers are making people easier to use everyday'
Examine what is said, not who speaks -- Silence betokens consent -- Love the truth but pardon error.
"Science is about questioning the status quo. Questioning authority".
In the absence of evidence, opinion is indistinguishable from prejudice.
Thanks a lot! that was what i was looking for
Re^2: Help with my Code
by madM (Beadle) on Aug 10, 2013 at 13:29 UTC
Thats what im trying to solve.. i schould get 44.530219 but im getting also 44.51357658 and the mutationRate is completly different :(

Create A New User
Node Status?
node history
Node Type: note [id://1048901]
help
Chatterbox?
and all is quiet...

How do I use this? | Other CB clients
Other Users?
Others pondering the Monastery: (4)
As of 2017-08-17 00:33 GMT
Sections?
Information?
Find Nodes?
Leftovers?
Voting Booth?
Who is your favorite scientist and why?

Results (277 votes). Check out past polls.

Notices?