help finding mean, numbers above and below mean from an exponential distribution

kaupifalco has asked for the wisdom of the Perl Monks concerning the following question:

Dear Monjes!!

I'm trying to generate a couple of things here. First I wrote a program that allows me to determine the generations that is needed for two lineages to meet back in time (evolving lineages). I'm assuming a mean of 4N/(k(k-1), with k= 2 and N=10000. My program is working and this is the script:

$k = 25; 
$N = 10000;

while ($k>=2)
    {
    if ($k==1)     
        {last;}
        $mean = (4*$N)/($k*($k-1));
        $time = (-log(rand)*$mean);
        
        push (@cls, $time);
        
        $k=$k-1;
        }
    open (GENE,'>coal.xls'); 
    for ($time=0; $time<24; $time ++) {
        print GENE "@cls[$time]\n";
        }
[download]

Now I need to write a simulation using an exponential distribution built from 5000 numbers, and afterwards get the mean, and the number of estimates above and below the mean. For this I need to use the same program as above. Although it is not working. So far this is my non-working scrypt:

#!usr/bin/perl
use warnings;

for ($sample=0; $sample<5000; $sample ++)
    {
    $k = 25; 
    $N = 10000;
    @cls=();
    $total=0;
    
while ($k>=2)
    {
    if ($k==1)     
        {last;}
        $mean = (4*$N)/($k*($k-1));
        $time = (-log(rand)*$mean);
        
        push (@cls, $time);
        $k=$k-1;
    }
        
        $ntotal+=$_ for @first;
        $nmean = $ntotal/5000;
        print "$nmean\n";
        foreach $first(@first)
        {
        if ($first<=$nmean) {push(@low, $first)}}
            $low= scalar @low;
            $hig= scalar (5000-@low);
                
            print "$low\n";
            print "$hig\n";
                
                
                }
[download]

Thanks in advance!

Comment on help finding mean, numbers above and below mean from an exponential distribution Select or Download Code

Replies are listed 'Best First'.
Re: help finding mean, numbers above and below mean from an exponential distribution by almut (Canon) on Oct 13, 2009 at 02:18 UTC
You're never assigning anything to `@first`, so it's no big surprise that the code using it doesn't do what you want... Also, do yourself a favor and use some consistent indenting; and of course strictures.	[reply] [d/l]
Re: help finding mean, numbers above and below mean from an exponential distribution by GrandFather (Saint) on Oct 13, 2009 at 04:07 UTC
By adding use strict, fixing the errors flagged by strictures, fixing the indentation and using Perl rather than C I get the following (which may or may not do what you want): `#!usr/bin/perl use warnings; use strict; my $ntotal; for (1 .. 5000) { my $k = 25; my $N = 10000; my @cls; my $total = 0; while ($k >= 2) { my $mean = (4 * $N) / ($k * ($k - 1)); my $time = (-log (rand) * $mean); $ntotal += $time; push @cls, $time; --$k; } my $nmean = $ntotal / 5000; print "$nmean\n"; my @low = grep {$_ <= $nmean} @cls; my $nLow = scalar @low; my $hig = scalar (5000 - $nLow); print "$nLow\n"; print "$hig\n"; }` [download] Prints: `2.37332066086488 0 5000 6.21051651542152 1 4999 8.59464183062829 0 5000 ...` [download] True laziness is hard work	[reply] [d/l] [select]
Re: help finding mean, numbers above and below mean from an exponential distribution by biohisham (Priest) on Oct 13, 2009 at 11:11 UTC
In both your programs you have: `while ($k>=2) { if ($k==1) {last;}` [download] so the condition for entering the while loop is $k>=2 but inside the loop itself someone might guess that you're checking if $k==1 and it may not become clear that there's a decrement towards the while block end, for readability and logic flow it would be better if this is written down towards the end of the while block after `$k=$k-1;` for example, but again, since you don't want $k==1 to be considered the loop would exit automatically if you only write `while($k>2)` without the need to check if $k has reached to 1. compare: `while($k>2){ $k=$k-1; print "$k\n"; } #and while($k>=2){ $k=$k-1; #last if $k==1; print "$k\n"; }` [download] $total=0; is used only once, what is it intended to contain?, where @first is coming from in `$ntotal+=$_ for @first;` Another major point, using strictures can save you a lot of nasty repercussions arising from undiscovered bugs and for readability you may want to pick up an indentation that is intuitive or use Perltidy, here is an initial tutorial. Excellence is an Endeavor of Persistence. Chance Favors a Prepared Mind.	[reply] [d/l] [select]
A reply falls below the community's threshold of quality. You may see it by logging in.