#!/usr/bin/perl
# A script by Erika Edwards and Casey Dunn to jack knife outgroups

use strict;
srand;

#Initialize variables, check arguments, and parse arguments
if ( @ARGV < 3 ){
  &intro;
  exit;
}

my $nout = shift @ARGV;      #The number of outgroups
my $nreps = shift @ARGV;     #The number of jackknife replicates to perform
my $nameblock = shift @ARGV; #Name of the file with the PAUP block to put between each jackknife replicate
my $nslice = shift @ARGV;    #Optional preset of the number of outgroups to include in each replicate

if ( $nslice > $nout ) {
    print "The number of outgroups to include exceeded the indicated number ";
    print "of outgroups. Quitting.\n";
}

my $nameout = "hotgroup.nex";

#Open the output nexus file
open (OUT, ">$nameout") or die "can't open $nameout";

#Write the file header
print OUT "#nexus\n[!This file was generated by hotgroups with the following settings:\n number of outgroups: $nout\n number of reps: $nreps\n paup block: $nameblock]\n\n";

#Loop through the number of specified reps
for (my $rep=0; $rep < $nreps; $rep++){

    #Make an array of consecutive integers numbered 1 through the number of outgroups
    my @seed;
    for (my $i=1; $i <= $nout; $i++ ) {
       push (@seed, $i);
    }
    #Randomize the order of the elements in the array
    my @rand = randomize ( @seed );

    print OUT "begin paup\;\ndelete ";

    
    # Decide how many outgroups to include
    my $amount;
    if ($nslice > 0) {
        #If specified, use the number indicated by the user
        $amount = $nslice;
    }
    else{
        # Pick a random number of out groups to include that ranges between 1 and the number of outgroups -1
        $amount = randint ( $nout-1 );
    }

    # Loop through and print the specified number of random outgroup taxa
    for (my $j=0; $j<$amount; $j++){
       my $element = shift @rand;
       print OUT "$element ";
    }

    print OUT " /only\;\nend\;\n\n";

    #Open the file with the block to insert to execute each run
    open (BLOCK, "$nameblock") or die "can't open $nameblock";
    my $nice = $rep + 1;
    while (my $line = <BLOCK> ) {
      #Replace the magic tree name with consecutive names if the user 
      #wants the trees from each run in a different file
      $line =~ s/REPLACE.tre/tree$nice.tre/; 

      #Print the block line to the output file
      print OUT "$line";
    }
    
    print OUT "\n\n";


}

#Let the user know that something happened
print "Done.\n";


sub randomize {
#Take an array and randomize the order of its elements
    my @start;
    @start = @_;
   
    my @out;
    
    while ( @start ) {
        my $n = @start;
        my $r = &randint( $n ) - 1;
        my $element = splice (@start, $r, 1);
        push (@out, $element);
    }


    return @out;
}


sub randint {
#Calculate a random integer from 1 to the specified number
   my $upper= shift @_;
   my $size = 1/$upper;

   my $n;
   my $f = rand;
   while ($f > 0) {
      $n++;
      $f = $f - $size;
   }

   return $n;

}


sub intro {
     print "\nhotgroup\nA script by Erika Edwards and Casey Dunn to jackknife outgroups.\n";
     print "October 2004, Ecology and Evolutionary Biology, Yale University\n\n";
     print "Usage:\nhotgroup NUMOUTGROUPS NUMREPS BLOCK [NUMRETAIN]\n\n";
     print "Where:\n";
     print "NUMOUTGROUPS is the number of outgroups in your dataset\n";
     print "NUMREPS is the number of jackknife replicates to perform\n";
     print "BLOCK is the name of the file containing the PAUP block that hotgroup will\n";
     print "insert between replicates\n";
     print "NUMRETAIN is an optional argument that can be used to fix the number of\n";
     print "outgroups included in each replicate\n\n";  
}