view release on metacpan or  search on metacpan

bin/alnpi  view on Meta::CPAN

  ## COMPUTE PROBABILITY OF ALLELIC CONFIGURATION USING KARLIN MACGREGOR 1972
  $Num_alleles = scalar keys %alleles;
  $H = 1;
  my %a;
  map {$H -= (($_/$nseqs) ** 2)} values %alleles;
  if (!$slide_window and !$pairwise and $Theta_w > 0) {
    $Exp_num_alleles = 1 / $Theta_w;
    for my $i (1 .. ($#seqs )) {
      $Exp_num_alleles += 1 / ($Theta_w + $i);
    }
    $Exp_num_alleles *= $Theta_w; ## Ewen's Sampling formula
    map {$a{$_}++} values %alleles; ## THIS MAKES THE a_i
    my $log_prob   = $Num_alleles * log ($Theta_w);
    #      $log_prob  += Math::Gsl::Sf::sf_lnfact($Num_alleles);
    $log_prob  += Ramanujan_logfact($Num_alleles);
    
    foreach my $i  (0 .. ($nseqs - 1)) {
      $log_prob -= log ($Theta_w + $i);
    }
    foreach (keys %a) {
      #   $log_prob -= ($a{$_} * log ($_)) + Math::Gsl::Sf::sf_lnfact($a{$_});