branch-1-6/Bio/PhyloNetwork/RandomFactory.pm

   1 # $Id$
   2 #
   3 # Module for Bio::PhyloNetwork::RandomFactory
   4 #
   5 # Please direct questions and support issues to <bioperl-l@bioperl.org>
   6 #
   7 # Cared for by Gabriel Cardona <gabriel(dot)cardona(at)uib(dot)es>
   8 #
   9 # Copyright Gabriel Cardona
  10 #
  11 # You may distribute this module under the same terms as perl itself
  12
  13 # POD documentation - main docs before the code
  14
  15 =head1 NAME
  16
  17 Bio::PhyloNetwork::RandomFactory - Module to generate random
  18 Phylogenetic Networks
  19
  20 =head1 SYNOPSIS
  21
  22  use strict;
  23  use warnings;
  24
  25  use Bio::PhyloNetwork;
  26  use Bio::PhyloNetwork::RandomFactory;
  27
  28  # Will generate at random all the 66 binary tree-child phylogenetic
  29  # networks with 3 leaves
  30
  31  my $factory=Bio::PhyloNetwork::RandomFactory->new(-numleaves=>3,-norepeat=>1);
  32
  33  my @nets;
  34
  35  for (my $i=0; $i<66; $i++) {
  36    my $net=$factory->next_network();
  37    push @nets,$net;
  38    print "".(scalar @nets).": ".$net->eNewick()."\n";
  39  }
  40
  41 =head1 DESCRIPTION
  42
  43 Builds a random (binary tree-child) phylogenetic network each time
  44 next_network is called.
  45
  46 =head1 AUTHOR
  47
  48 Gabriel Cardona, gabriel(dot)cardona(at)uib(dot)es
  49
  50 =head1 SEE ALSO
  51
  52 L<Bio::PhyloNetwork>
  53
  54 =head1 APPENDIX
  55
  56 The rest of the documentation details each of the object methods.
  57
  58 =cut
  59
  60 package Bio::PhyloNetwork::RandomFactory;
  61
  62 use strict;
  63 use warnings;
  64
  65 use base qw(Bio::Root::Root);
  66
  67 use Bio::PhyloNetwork;
  68 use Math::Random;
  69 use Bio::Tree::RandomFactory;
  70
  71 =head2 new
  72
  73  Title   : new
  74  Usage   : my $factory = new Bio::PhyloNetwork::RandomFactory();
  75  Function: Creates a new Bio::PhyloNetwork::RandomFactory
  76  Returns : Bio::PhyloNetwork::RandomFactory
  77  Args    : -numleaves => integer
  78             OR
  79            -leaves => reference to an array (of leaves names)
  80            -numhybrids => integer [optional]
  81            -norepeat => boolean [optional]
  82
  83 Returns a Bio::PhyloNetwork::RandomFactory object. Such an object will create
  84 random binary tree-child phylogenetic networks each time next_network
  85 is called.
  86
  87 If the parameter -leaves=E<gt>\@leaves is given, then the set of leaves of
  88 these networks will be @leaves. If it is given the parameter
  89 -numleaves=E<gt>$numleaves, then the set of leaves will be "l1"..."l$numleaves".
  90
  91 If the parameter -numhybrids=E<gt>$numhybrids is given, then the generated
  92 networks will have exactly $numhybrids hybrid nodes. Note that, necessarily,
  93 $numhybrids E<lt> $numleaves. Otherwise, the number of hybrid nodes will be chosen
  94 at random for each call of next_network.
  95
  96 If the parameter -norepeat=E<gt>1 is given, then successive calls of next_network
  97 will give non-isomorphic networks.
  98
  99 =cut
 100
 101 sub new {
 102   my ($pkg,@args)=@_;
 103
 104   my $self=$pkg->SUPER::new(@args);
 105
 106   my ($leavesR,$numleaves,$numhybrids,$norepeat)=
 107     $self->_rearrange([qw(LEAVES
 108                           NUMLEAVES
 109                           NUMHYBRIDS
 110                           NOREPEAT)],@args);
 111   my @leaves;
 112   if ((! defined $leavesR) && (defined $numleaves)) {
 113     @leaves=map {"l$_"} (1..$numleaves);
 114     $leavesR=\@leaves;
 115   }
 116   if (! defined $leavesR) {
 117     $self->throw("No leaves set neither numleaves given");
 118   }
 119   $norepeat ||= 0;
 120
 121   $self->{leaves}=\@leaves;
 122   $self->{numleaves}=$numleaves;
 123   $self->{numhybrids}=$numhybrids if defined $numhybrids;
 124   $self->{norepeat}=$norepeat;
 125   $self->{found}=[];
 126   $self->{tree_factory}=Bio::Tree::RandomFactory->new(-taxa => \@leaves);
 127   bless($self,$pkg);
 128 }
 129
 130 =head2 next_network
 131
 132  Title   : next_network
 133  Usage   : my $net=$factory->next_network()
 134  Function: returns a random network
 135  Returns : Bio::PhyloNetwork
 136  Args    : none
 137
 138 =cut
 139
 140 sub next_network {
 141   my ($self)=@_;
 142
 143   my $numleaves=$self->{numleaves};
 144   my @found=@{$self->{found}};
 145   my $numhybrids;
 146  START:
 147   if (! defined $self->{numhybrids}) {
 148     $numhybrids=int(rand($numleaves));
 149   }
 150   else {
 151     $numhybrids=$self->{numhybrids};
 152   }
 153   my $tf=$self->{tree_factory};
 154   my $tree=$tf->next_tree;
 155   my $net=Bio::PhyloNetwork->new(-tree=>$tree);
 156   for (my $i=1; $i<=$numhybrids; $i++) {
 157     $net=random_attack($net,$i);
 158   }
 159   if ($self->{norepeat}) {
 160     foreach my $ant (@found) {
 161       goto START if $net->is_mu_isomorphic($ant);
 162     }
 163     push @found,$net;
 164     $self->{found}=\@found;
 165   }
 166   return $net;
 167 }
 168
 169 sub random_attack {
 170   my ($net,$lbl)=@_;
 171
 172   my $graph=$net->{graph};
 173   my ($u1,$v1,$u2,$v2);
 174   do {
 175     my $e1=$graph->random_edge;
 176     my $e2=$graph->random_edge;
 177     $u1=$e1->[0];
 178     $v1=$e1->[1];
 179     $u2=$e2->[0];
 180     $v2=$e2->[1];
 181   } while (! $net->is_attackable($u1,$v1,$u2,$v2,$lbl));
 182   $net->do_attack($u1,$v1,$u2,$v2,$lbl);
 183   return $net;
 184 }
 185
 186 1;