#!/usr/local/bin/perl
# GA.pl
# defaults given are those in the paper
=head1 NOTES
As with all the scripts, Perl must be installed.
For this script alone, (but none of the others) the Tk library must also be in the lib.
On Linux (or similar) you'll need to get Tk from CPAN (www.cpan.org) and install it yourself.
On Windows, the best option is to get a recent version of ActivePerl (www.activestate.com), which
is free, and has all the Tk libraries ready to go.
The usage is:
perl -w GA.pl
As with all Perl scripts, they are offered as is, without any warranties concerning
safety, functionality etc
The author welcomes any comments on any aspects of the software.
to: dgatherer2002@yahoo.co.uk
This is an occasionally on-going project and newer versions of code are frequently
available, so if you want the really rough cuts, please ask.
=cut
#------PRELIMINARIES-------------------------------------------------#
use strict;
use POSIX;
srand(); # initialise randomiser
# user-defined global variables are:
# 1. length of string (all strings same length in this version)
# 2. size of population
# 3. number of generations to run
# 4. mutation rate
# 5. recombination rate
# 6. lower cut-off for negative selection (in st.devs.)
# 7. upper cut-off for mutation (to protect best strings from mut)
# 8. upper cut-off for recombination (likewise from recomb.)
# 1-3 must be positive integers
# 4 and 5 must be positive floats
# 6-8 are floats
# all have defaults, as follows:
# DEFAULTS FOR GA
my $GA_len = 6; # ie. reprate, migrate, contage, nn_cont, cultsel, ysel
my $GA_popsize = 50;
my $GA_time = 15; # 50 strings of 6, over 15 iterations
my $GA_mutrat = 0.05; # low mutation
my $GA_recrat = 0.1; # recombine (subject to no protection)
my $discard = 0.775; # more than part of a standard dev. below average? then die
my $no_mut = 0.5; # where theshold is applied for mut and rec
my $no_rec = 0.5; # all below average mutate or recombine
my $attrib; # used in recomb etc.
my $direction = "D"; # U for up ie. max isolation selected for
my @the_G_string = qw(rep_rate
contag_rate
nn_contag_rate
mig_rate
cult_sel
nat_sel); # must be 'strung' although initialised as hash
# GLOBALS FOR SELECTION
my $popsize = 200; # initial pop of 200
my $time = 100; # over 100 gens
my $line = my $cols = 10; # in a 10 by 10 landscape
my @pop;
my $pause = 1; # when to stop and draw the map
my $capacity = 100000000;
my $reps = 2;
open GRAPH, ">GA_graph.xls" or die "\ncan't open the graph";
print GRAPH "\nstrings: $GA_popsize over $GA_time";
print GRAPH "\nmutrat: $GA_mutrat";
print GRAPH "\nrecrat: $GA_recrat";
print GRAPH "\nthresholds: die: $discard, recom: $no_rec, mut: $no_mut";
print GRAPH "\nAssessment time: $time";
print GRAPH "\nGen.\tAv.fitness\tst.dev.";
#------MAIN CODE SEQUENCE---------------------------------------------#
my $GA_pop = &GA_initialise(); # make a population of strings
print "\ninitialised...";
&GA_pop_state($GA_pop); # look at it
for(my $gen=1;$gen<=$GA_time;$gen++) # for a certain number of cycles
{
print "\nGeneration $gen";
print GRAPH "\n$gen";
$GA_pop = &GA_assess($GA_pop); # perform simulation calculation
$GA_pop = &GA_select($GA_pop); # perform fitness assessment
print "\nselected...";
&GA_pop_state($GA_pop);
$GA_pop = &GA_reproduce($GA_pop); # discard bad, replicate good
# &GA_pop_state($GA_pop); # look at it
$GA_pop = &GA_recombine($GA_pop); # recombine it
print "\nrecombined...";
# &GA_pop_state($GA_pop);
$GA_pop = &GA_mutate($GA_pop); # mutate it
print "\nmutated...";
# $GA_pop = &GA_select($GA_pop);
# &GA_pop_state($GA_pop);
}
print "\nresults of GA process";
$GA_pop = &GA_assess($GA_pop); # perform simulation calculation
$GA_pop = &GA_select($GA_pop);
&GA_pop_state($GA_pop);
#-----SUBROUTINES-----------------------------------------------------#
sub GA_initialise() # initialise population
{
my @GA_array; # initialise array for GA
my $av_fitness = 0; # initially zero fitness
my $st_dev = 0;
for(my $y=0;$y<=$GA_popsize-1;$y++) # run through GA pop
{
my %one_string; # the GA string is actually a hash
$one_string{rep_rate} = 0.01;
$one_string{contag_rate} = rand();
$one_string{nn_contag_rate} = rand();
$one_string{mig_rate} = 0.1*rand(); # don't want ultra-high mig
$one_string{cult_sel} = rand();
$one_string{nat_sel} = rand();
$one_string{g_disp} = 0;
$one_string{t_disp} = 0;
$one_string{fitness} = $one_string{t_disp};
($one_string{cult_sel} >=0.5) ? ($one_string{cult_sel} = "Y") : ($one_string{cult_sel} = "N");
($one_string{nat_sel} >=0.5) ? ($one_string{nat_sel} = "Y") : ($one_string{nat_sel} = "N");
push @GA_array, \%one_string;
}
# RECORD OF POP
my $new_pop =
{
POPULATION => [ @GA_array ],
FITNESS => $av_fitness,
STDEV => $st_dev,
};
return $new_pop;
}
#--------------------------------------------------------------------
sub GA_recombine() # rather similar to mutate
{
print "\nin recombine\n";
my $input = $_[0];
my $av_fitness = $input ->{FITNESS};
my $st_dev = $input ->{STDEV};
my $threshold = $av_fitness - ($no_rec * $st_dev); # threshold at $discard st_devs above mean
my @current_pop = @{ $input ->{POPULATION} };
>
my @recombed_pop;
for(my $i=0; $i<=$#current_pop; $i++)
{
my %indiv = %{$current_pop[$i]};
my $rand = rand();
if(($rand <= $GA_recrat)) # && ($indiv{t_disp} <= $threshold)) # only recombine subject to recombination rate (default always)
{
my $chosen_one = -1; # initialise out of range
while($chosen_one < 0)
{
my $rand = rand(); # can recombine with self
$chosen_one = floor($rand*$GA_popsize);
}
print "\nchosen: ".($chosen_one+1)." with: ".($i+1);
my %other = %{$current_pop[$chosen_one]};
my $rand = rand();
my $synapse = floor($rand*$GA_len);
print " at ".($synapse+1);
my %new_GA;
for(my $a=0;$a<=$synapse;$a++)
{
my $attrib = $the_G_string[$a];
$new_GA{$attrib} = $indiv{$attrib};
}
for(my $a=$synapse+1;$a<=$GA_len-1;$a++)
{
my $attrib = $the_G_string[$a];
$new_GA{$attrib} = $other{$attrib};
}
$new_GA{t_disp} = $other{t_disp};
$new_GA{g_disp} = $other{g_disp};
$new_GA{fitness} = $other{fitness};
push @recombed_pop, \%new_GA;
}
else { push @recombed_pop, \%indiv; }
}
# RECORD OF POP
my $new_pop =
{
POPULATION => [ @recombed_pop ],
FITNESS => $av_fitness,
STDEV => $st_dev,
};
return $new_pop;
}
#----------------------------------------------------------------------
sub GA_pop_state() # prints current pop
{
my $input = $_[0];
my $av_fitness = $input ->{FITNESS};
my $st_dev = $input ->{STDEV};
$st_dev = sprintf("%.2f", $st_dev);
my @current_pop = @{ $input ->{POPULATION} };
>
print "\nfitness: $av_fitness stdev: $st_dev\n";
print GRAPH "\t$av_fitness\t$st_dev";
for(my $p=0; $p<=$#the_G_string; $p++)
{
print "$the_G_string[$p] ";
}
print "t_disp g_disp f";
my $line = 0;
foreach(@current_pop)
{
my %indiv = %$_;
print "\n".++$line." ";
for(my $p=0; $p<=$#the_G_string-2; $p++)
{
my $attrib = $the_G_string[$p];
printf("%.2f\t", $indiv{$attrib});
}
for(my $p=$#the_G_string-1; $p<=$#the_G_string; $p++)
{
my $attrib = $the_G_string[$p];
print "\t$indiv{$attrib}";
}
if(exists($indiv{t_disp}) && exists($indiv{g_disp}))
{
print "\t$indiv{t_disp}\t$indiv{g_disp}"; # \t$indiv{fitness}";
}
}
}
#------------------------------------------------------------------------
sub GA_mutate() # rather similar to initialise
{
print "\nin GA_mutate\n";
my $input = $_[0];
my $av_fitness = $input ->{FITNESS};
my $st_dev = $input ->{STDEV};
my $threshold = $av_fitness - ($no_mut*$st_dev); # threshold at 2 st_devs below mean (lowest 5%)
my @muted_pop;
my @current_pop = @{ $input ->{POPULATION} };
>
for(my $i=0; $i<=$#current_pop; $i++)
{
my %indiv = %{$current_pop[$i]};
my $rand = rand();
if(($rand <= $GA_mutrat)) ## && ($indiv{t_disp} <= $threshold)) # only mutate subject to rate (default always)
{
$rand = rand();
my $mut_point = floor($rand*$GA_len);
print "\nat ".($mut_point+1);
print "in $i going to mutate $the_G_string[$mut_point]";
my %new_GA;
my $mutated = $the_G_string[$mut_point];
for(my $p=0; $p<=$#the_G_string; $p++)
{
$attrib = $the_G_string[$p];
$new_GA{$attrib} = $indiv{$attrib};
}
if($mut_point <= 3)
{
$new_GA{$mutated} = $rand;
}
elsif($mut_point == 4)
{
($rand >=0.5) ? ($new_GA{cult_sel} = "Y") : ($new_GA{cult_sel} = "N");
}
elsif($mut_point == 5)
{
($rand >=0.5) ? ($new_GA{nat_sel} = "Y") : ($new_GA{nat_sel} = "N");
}
$new_GA{t_disp} = $indiv{t_disp};
$new_GA{g_disp} = $indiv{g_disp};
$new_GA{fitness} = $indiv{fitness};
$new_GA{rep_rate} = 0.01; # can't mutate rep rate
# $new_GA{mig_rate} = 0.01; # can't mutate mig rate, restor if you want to
push @muted_pop, \%new_GA;
}
else { push @muted_pop, \%indiv; }
}
# RECORD OF POP
my $new_pop =
{
POPULATION => [ @muted_pop ],
FITNESS => $av_fitness,
STDEV => $st_dev,
};
return $new_pop;
}
#------------------------------------------------------------------------
sub GA_select()
{
print "\nin GA_select\n";
my $total_fitness = 0; # average population fitness
my $input = $_[0];
my @current_pop = @{ $input ->{POPULATION} };
>
foreach(@current_pop)
{
my %indiv = %$_;
$total_fitness += $indiv{fitness};
}
my $av_fitness = $total_fitness/$GA_popsize;
my $fit_variance = 0; # variance of fitness
my $deviation = 0; # deviation from mean
my $st_dev = 0; # standard deviation
foreach(@current_pop)
{
my %indiv = %$_;
$deviation += ($indiv{fitness} - $av_fitness)**2;
}
$fit_variance = $deviation/($GA_popsize-1);
$st_dev = sqrt($fit_variance);
# RECORD OF POP
my $new_pop =
{
POPULATION => [ @current_pop ],
FITNESS => $av_fitness,
STDEV => $st_dev,
};
return $new_pop;
}
#--------------------------------------------------------------------
sub GA_reproduce() # rather similar to recombine
{
print "\nin reproduce\n";
my $input = $_[0];
my $av_fitness = $input ->{FITNESS};
my $st_dev = $input ->{STDEV};
my $threshold; # threshold at 2 st_devs below mean (lowest 5%)
my @current_pop = @{ $input ->{POPULATION} };
>
my @selected_pop;
# discard all those under specified fitness
for(my $i=0; $i<=$#current_pop; $i++)
{
my %indiv = %{$current_pop[$i]};
if($direction eq "U") # ie seln for max isol
{
$threshold = ($av_fitness - ($discard*$st_dev));
if($indiv{t_disp} <= $threshold) # only discard if score is under average
{
my $chosen_one = -1; # initialise out of range
while($chosen_one < 0)
{
my $rand = rand(); # can recombine with self
$chosen_one = floor($rand*$GA_popsize);
}
print "\nchosen: ".($chosen_one+1)." replacing: ".($i+1);
my %other = %{$current_pop[$chosen_one]};
my %new_GA;
for(my $p=0; $p<=$#the_G_string; $p++)
{
$attrib = $the_G_string[$p];
$new_GA{$attrib} = $other{$attrib};
}
$new_GA{t_disp} = $other{t_disp};
$new_GA{g_disp} = $other{g_disp};
$new_GA{fitness} = $other{g_disp};
push @selected_pop, \%new_GA;
}
else { push @selected_pop, \%indiv; } # keep it
}
elsif($direction eq "D") # ie seln for min isol
{
$threshold = ($av_fitness + ($discard*$st_dev));
if($indiv{t_disp} >= $threshold) # only discard if score is under average
{
my $chosen_one = -1; # initialise out of range
while($chosen_one < 0)
{
my $rand = rand(); # can recombine with self
$chosen_one = floor($rand*$GA_popsize);
}
print "\nchosen: ".($chosen_one+1)." replacing: ".($i+1);
my %other = %{$current_pop[$chosen_one]};
my %new_GA;
for(my $p=0; $p<=$#the_G_string; $p++)
{
$attrib = $the_G_string[$p];
$new_GA{$attrib} = $other{$attrib};
}
$new_GA{t_disp} = $other{t_disp};
$new_GA{g_disp} = $other{g_disp};
$new_GA{fitness} = $other{g_disp};
push @selected_pop, \%new_GA;
}
else { push @selected_pop, \%indiv; } # keep it
}
}
# RECORD OF POP
my $new_pop =
{
POPULATION => [ @selected_pop ],
FITNESS => $av_fitness,
STDEV => $st_dev,
};
return $new_pop;
}
#----SUBROUTINES-----------------------------------------------#
sub GA_assess
{
print "\nin GA_assess\n";
my $input = $_[0];
my $av_fitness = $input ->{FITNESS};
my $st_dev = $input ->{STDEV};
my @current_pop = @{ $input ->{POPULATION} };
>
my @new_pop;
foreach(@current_pop)
{
my %indiv = %$_;
for(my $x=0; $x<=$line-1; $x++)
{
for(my $y=0; $y<=$cols-1; $y++)
{
$pop[$x][$y] = 0;
}
}
my $g_disp_av = my $t_disp_av = 0;
for(my $runs=1; $runs<=$reps; $runs++)
{
$popsize = 200; # reset for each run
my $pop = &evenly_initialise();
my $g_disp = my $t_disp = 0;
for(my $gen=1;$gen<=$time;$gen++) # for a certain number of cycles
{

$pop = &reproduce($pop, \%indiv); # replicate and pass traits by contagion
$pop = &migrate($pop, \%indiv); # move strings over grid
if($gen%$pause ==0)
{
($g_disp, $t_disp) = &pop_state($pop);
} # look at it
}
$g_disp_av += $g_disp;
$t_disp_av += $t_disp;
}
$g_disp_av /= $reps;
$t_disp_av /= $reps;
$indiv{g_disp} = $g_disp_av;
$indiv{t_disp} = $t_disp_av;
$indiv{fitness} = $t_disp_av;
push @new_pop, \%indiv;
}
my $new_pop =
{
POPULATION => [ @new_pop ],
FITNESS => $av_fitness,
STDEV => $st_dev,
};
return $new_pop;
}
#----------------------------------------------------------------------
sub pop_state() # prints current pop
{
my $input = $_[0];
my $av_fitness = $input ->{FITNESS};
my @current_pop = @{ $input ->{POPULATION} };
>
my @gene_state; # arrays for calculating the
my @trait_state; # horizontal and vertical dispersal
for(my $x=0; $x<=$line-1; $x++)
{
for(my $y=0; $y<=$cols-1; $y++)
{
my $freq_A = my $freq_1 = 0;
my $cell_count = 0; # how many per cell
foreach(@current_pop)
{
my $xcoord = $_ -> { LOCATION }[0];
my $ycoord = $_ -> { LOCATION }[1];
if($xcoord == $x && $ycoord == $y)
{
$cell_count++;
if($_ -> { TRAIT } eq "A")
{
$freq_A++;
}
if($_ -> { SERIAL_NUM } == 1)
{
$freq_1++;
}
}
}
$pop[$x][$y] = $cell_count;
if($freq_A == 0 && $freq_1 == 0 && $pop[$x][$y] == 0)
{
$trait_state[$x][$y] = 0;
$gene_state[$x][$y] = 0;
}
else
{
$freq_A /= $pop[$x][$y];
$freq_1 /= $pop[$x][$y];
$trait_state[$x][$y] = $freq_A;
$gene_state[$x][$y] = $freq_1;
}
}
}
my($g, $t) = &disp_calc(\@gene_state, \@trait_state);
return ($g, $t);
}
#--------------------------------------------------------------------
sub reproduce() # rather similar to recombine
{
my ($input, $params) = @_;
my %params = %$params;
my $rep_rate = $params{rep_rate};
my $contag_rate = $params{contag_rate};
my $nn_contag_rate = $params{nn_contag_rate};
my $cult_sel = $params{cult_sel};
my $nat_sel = $params{nat_sel};
my $contag_events = my $repro_events = 0; # how many per cycles
for(my $b=0;$b<=$popsize-1;$b++) # run through pop
{
my $rand = rand();

$rand = rand();
if($nat_sel =~ /Y/i && @{ $input ->{POPULATION} }[$b] ->{TRAIT} eq "A")
{
$rand /= 2; # double the random number
}
if($rand < $rep_rate) # also duplicate every 10th one
{
my $xcoord = @{ $input ->{POPULATION} }[$b] -> { LOCATION }[0];
my $ycoord = @{ $input ->{POPULATION} }[$b] -> { LOCATION }[1];
if($pop[$xcoord][$ycoord] < $capacity)
{
my @new_one_seq;
my $new_seq =
{
LOCATION => [ $xcoord, $ycoord ],
FITNESS => @{ $input ->{POPULATION} }[$b] ->{FITNESS},
SERIAL_NUM => @{ $input ->{POPULATION} }[$b] ->{SERIAL_NUM},
LINEAGE => @{ $input ->{POPULATION} }[$b] ->{LINEAGE},
TRAIT => @{ $input ->{POPULATION} }[$b] ->{TRAIT},
};
@{ $input ->{POPULATION} }[$popsize++] = $new_seq;
$pop[$xcoord][$ycoord]++;
$repro_events++;
}
}
$rand = rand();
if($cult_sel =~ /Y/i && @{ $input ->{POPULATION} }[$b] ->{TRAIT} eq "A")
{
$rand *= 2; # halve the random number
}
if($rand < $contag_rate) # and contage every 10th one

{
my $chosen_one = floor($rand*$popsize);
while(abs(@{ $input ->{POPULATION} }[$b] -> { LOCATION }[0]
- @{ $input ->{POPULATION} }[$chosen_one] ->{ LOCATION }[0]) > 1
|| abs(@{ $input ->{POPULATION} }[$b] -> { LOCATION }[1]
- @{ $input ->{POPULATION} }[$chosen_one] -> { LOCATION }[1]) > 1)
{ # replace with random other string
my $chos_rand = rand();
$chosen_one = floor($chos_rand*$popsize);
}
# first allow (smaller?) chance of contagion to neighbours
if($rand < $nn_contag_rate)
{
if(abs(@{ $input ->{POPULATION} }[$b] -> { LOCATION }[0]
- @{ $input ->{POPULATION} }[$chosen_one] ->{ LOCATION }[0]) <= 1
&& abs(@{ $input ->{POPULATION} }[$b] -> { LOCATION }[1]
- @{ $input ->{POPULATION} }[$chosen_one] -> { LOCATION }[1]) <= 1)
{
@{ $input ->{POPULATION} }[$b] ->{TRAIT} = @{ $input ->{POPULATION}} }[$chosen_one] ->{TRAIT};
my $xcoord = @{ $input ->{POPULATION} }[$b] -> { LOCATION }[0]; # direct access, better ??
my $ycoord = @{ $input ->{POPULATION} }[$b] -> { LOCATION }[1]; # yes, same method for string access
my $trait = @{ $input ->{POPULATION} }[$b] -> { TRAIT };
$contag_events++;
next; # only one contag per indiv per gen
}
}
else
{
while(abs(@{ $input ->{POPULATION} }[$b] -> { LOCATION }[0]
- @{ $input ->{POPULATION} }[$chosen_one] ->{ LOCATION }[0]) > 0
|| abs(@{ $input ->{POPULATION} }[$b] -> { LOCATION }[1]
- @{ $input ->{POPULATION} }[$chosen_one] -> { LOCATION }[1]) > 0)
{ # replace with random other string
my $chos_rand = rand();
$chosen_one = floor($chos_rand*$popsize);
}
@{ $input ->{POPULATION} }[$b] ->{TRAIT} = @{ $input ->{POPULATION}} }[$chosen_one] ->{TRAIT};
my $xcoord = @{ $input ->{POPULATION} }[$b] -> { LOCATION }[0]; # direct access, better ??
my $ycoord = @{ $input ->{POPULATION} }[$b] -> { LOCATION }[1]; # yes, same method for string access
my $trait = @{ $input ->{POPULATION} }[$b] -> { TRAIT };
$contag_events++;
next;
}
}
}
return $input;
}
#--------------------------------------------------------------------
sub migrate() # rather similar to reproduce
{
my ($input, $params) = @_;
my %params = %$params;
my $mig_rate = $params{mig_rate};
my $migrats = 0; # how many
for(my $b=0;$b<=$popsize-1;$b++)
{
my $rand = rand();
if($rand < $mig_rate)
{ # choose
$rand = rand(); # random direction
my $xcoord = @{ $input ->{POPULATION} }[$b] ->{LOCATION}[0];
my $ycoord = @{ $input ->{POPULATION} }[$b] ->{LOCATION}[1];
if($rand <= 0.125 && $xcoord >0 && $ycoord >0) # go up and left
{
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[0]--;
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[1]--;
$pop[$xcoord][$ycoord]--;
$pop[$xcoord-1][$ycoord-1]++;
$migrats++;
}
elsif($rand <= 0.25 && $xcoord >0) # go straight up
{
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[0]--;
$pop[$xcoord-1][$ycoord]++;
$pop[$xcoord][$ycoord]--;
$migrats++;
}
elsif($rand <= 0.375 && $xcoord >0 && $ycoord <$cols-1) # go up and right
{
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[0]--;
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[1]++;
$pop[$xcoord-1][$ycoord+1]++;
$pop[$xcoord][$ycoord]--;
$migrats++;
}
elsif($rand <= 0.5 && $ycoord <$cols-1) # go straight right
{
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[1]++;
$pop[$xcoord][$ycoord+1]++;
$pop[$xcoord][$ycoord]--;
$migrats++;
}
elsif($rand <= 0.625 && $xcoord <$line-1 && $ycoord <$cols-1) # go down and right
{
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[0]++;
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[1]++;
$pop[$xcoord+1][$ycoord+1]++;
$pop[$xcoord][$ycoord]--;
$migrats++;
}
elsif($rand <= 0.75 && $xcoord <$line-1) # go straight down
{
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[0]++;
$pop[$xcoord+1][$ycoord]++;
$pop[$xcoord][$ycoord]--;
$migrats++;
}
elsif($rand <= 0.825 && $xcoord < $line-1 && $ycoord > 0) # go down and left
{
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[0]++;
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[1]--;
$pop[$xcoord+1][$ycoord-1]++;
$pop[$xcoord][$ycoord]--;
$migrats++;
}
elsif($ycoord >0) # go straight left
{
@{ $input ->{POPULATION} }[$b] ->{LOCATION}[1]--;
$pop[$xcoord][$ycoord-1]++;
$pop[$xcoord][$ycoord]--;
$migrats++;
}
}
}
return $input;
}
#--------------------------------------------------------------------------------------------#
sub evenly_initialise() # initialise population
{
my @seq_array; # initialise array for seqs
my $indiv_fitness = 0; # fitness of string, only used when seln. operating
my $av_fitness = 0; # for whole pop, only used when seln. operating
my $st_dev = 0; # ditto
for(my $x=0; $x<=$line-1; $x++)
{
for(my $y=0; $y<=$cols-1; $y++)
{
if($popsize < ($line*$cols)) { die "popsize too small" };
for(my $z=0;$z<=($popsize/($line*$cols))-1;$z++) # run through pop
{
my @one_seq; # re-initialise each string
my $rand = rand(); # decide on horiz trait state
my $horiz_trait;
my $initial_id;
if ($rand < 0.25)
{
#print "A" ;
$horiz_trait = "A";
}
elsif ($rand < 0.5)
{
#print "B" ;
$horiz_trait = "B";
}
elsif ($rand < 0.75)
{
#print "C" ;
$horiz_trait = "C";
}
else
{
#print "O";
$horiz_trait = "O";
}
$rand = rand(); # decide on vertical trait state
if ($rand < 0.25)
{
#print "A" ;
$initial_id = 1;
}
elsif ($rand < 0.5)
{
#print "B" ;
$initial_id = 2;
}
elsif ($rand < 0.75)
{
#print "C" ;
$initial_id = 3;
}
else
{
#print "O";
$initial_id = 4;
}
# RECORD OF SEQUENCE
my $new_seq =
{
LOCATION => [ $x, $y ], # all are initially in square 0,0
FITNESS => $indiv_fitness, # if required
SERIAL_NUM => $initial_id, # for lineage mapping
LINEAGE => $initial_id, # will be grown, as lineage unfolds
TRAIT => $horiz_trait, # will pass by contagion
};
push @seq_array, $new_seq;
}
}
}
# RECORD OF POP
my $new_pop =
{
POPULATION => [ @seq_array ],
FITNESS => $av_fitness, # if required
STDEV => $st_dev, # if required
};
return $new_pop;
}
#--------------------------------------------------------------------------------------------#
sub disp_calc() # calculates dispersals
{
my($gene_state, $trait_state) = @_;
my $total_genes = my $total_traits = my $nbrd_genes = my $nbrd_traits =0;
my @genes = @$gene_state;
my @traits = @$trait_state;
for(my $x=0; $x<=$line-1; $x++)
{
for(my $y=0; $y<=$cols-1; $y++)
{
if($genes[$x][$y] >= 0.5)
{
$total_genes++;
unless($x == 0 || $y == 0 || $x == $line-1 || $y == $cols-1)
{
if($genes[$x+1][$y] >= 0.5
|| $genes[$x-1][$y] >= 0.5
|| $genes[$x][$y+1] >= 0.5

|| $genes[$x][$y-1] >= 0.5
|| $genes[$x+1][$y+1] >= 0.5
||$genes[$x-1][$y-1] >= 0.5
|| $genes[$x+1][$y-1] >= 0.5
|| $genes[$x-1][$y+1] >= 0.5)
{
$nbrd_genes++;
}
}
else
{
if($x == 0 && $y == 0) # top left
{
if($genes[$x+1][$y] >= 0.5
|| $genes[$x][$y+1] >= 0.5
|| $genes[$x+1][$y+1] >= 0.5)
{
# print "\tnn also maj.";
$nbrd_genes++;
}
}
elsif($x == 0 && $y == $cols-1) # top right
{
if($genes[$x+1][$y] >= 0.5
|| $genes[$x][$y-1] >= 0.5
|| $genes[$x+1][$y-1] >= 0.5)
{
# print "\tnn also maj.";
$nbrd_genes++;
}
}
elsif($x == $line-1 && $y == 0) # bottom left
{
if($genes[$x-1][$y] >= 0.5
|| $genes[$x][$y+1] >= 0.5
|| $genes[$x-1][$y+1] >= 0.5)
{
$nbrd_genes++;
}
}
elsif($x == $line-1 && $y == $cols-1) # bottom right
{
if($genes[$x-1][$y] >= 0.5
|| $genes[$x][$y-1] >= 0.5
|| $genes[$x-1][$y-1] >= 0.5)
{
$nbrd_genes++;
}
}
elsif($x == 0) # top row
{
if($genes[$x+1][$y] >= 0.5
|| $genes[$x][$y+1] >= 0.5
|| $genes[$x][$y-1] >= 0.5
|| $genes[$x+1][$y+1] >= 0.5
|| $genes[$x+1][$y-1] >= 0.5)
{
$nbrd_genes++;
}
}
elsif($y == 0) # left column
{
if($genes[$x+1][$y] >= 0.5
|| $genes[$x-1][$y] >= 0.5
|| $genes[$x][$y+1] >= 0.5
|| $genes[$x+1][$y+1] >= 0.5
|| $genes[$x-1][$y+1] >= 0.5)
{
$nbrd_genes++;
}
}
elsif($x == $line-1) # bottom row
{
if($genes[$x-1][$y] >= 0.5
|| $genes[$x][$y+1] >= 0.5
|| $genes[$x][$y-1] >= 0.5
||$genes[$x-1][$y-1] >= 0.5
|| $genes[$x-1][$y+1] >= 0.5)
{
$nbrd_genes++;
}
}
elsif($y == $cols-1) # right column
{
if($genes[$x+1][$y] >= 0.5
|| $genes[$x-1][$y] >= 0.5
|| $genes[$x][$y-1] >= 0.5
||$genes[$x-1][$y-1] >= 0.5
|| $genes[$x+1][$y-1] >= 0.5)
{
$nbrd_genes++;
}
}
}
}
if($traits[$x][$y] >= 0.5)
{
$total_traits++;
unless($x == 0 || $y == 0 || $x == $line-1 || $y == $cols-1)
{
if($traits[$x+1][$y] >= 0.5
|| $traits[$x-1][$y] >= 0.5
|| $traits[$x][$y+1] >= 0.5
|| $traits[$x][$y-1] >= 0.5
|| $traits[$x+1][$y+1] >= 0.5
||$traits[$x-1][$y-1] >= 0.5
|| $traits[$x+1][$y-1] >= 0.5
|| $traits[$x-1][$y+1] >= 0.5)
{
$nbrd_traits++;
}
}
else
{
if($x == 0 && $y == 0) # top left
{
if($traits[$x+1][$y] >= 0.5
|| $traits[$x][$y+1] >= 0.5
|| $traits[$x+1][$y+1] >= 0.5)
{
$nbrd_traits++;
}
}
elsif($x == 0 && $y == $cols-1) # top right
{
if($traits[$x+1][$y] >= 0.5
|| $traits[$x][$y-1] >= 0.5
|| $traits[$x+1][$y-1] >= 0.5)
{
$nbrd_traits++;
}
}
elsif($x == $line-1 && $y == 0) # bottom left
{
if($traits[$x-1][$y] >= 0.5
|| $traits[$x][$y+1] >= 0.5
|| $traits[$x-1][$y+1] >= 0.5)
{
$nbrd_traits++;
}
}
elsif($x == $line-1 && $y == $cols-1) # bottom right
{
if($traits[$x-1][$y] >= 0.5
|| $traits[$x][$y-1] >= 0.5
|| $traits[$x-1][$y-1] >= 0.5)
{
$nbrd_traits++;
}
}
elsif($x == 0) # top row
{
if($traits[$x+1][$y] >= 0.5
|| $traits[$x][$y+1] >= 0.5
|| $traits[$x][$y-1] >= 0.5
|| $traits[$x+1][$y+1] >= 0.5
|| $traits[$x+1][$y-1] >= 0.5)
{
$nbrd_traits++;
}
}
elsif($y == 0) # left column
{
if($traits[$x+1][$y] >= 0.5
|| $traits[$x-1][$y] >= 0.5
|| $traits[$x][$y+1] >= 0.5
|| $traits[$x+1][$y+1] >= 0.5
|| $traits[$x-1][$y+1] >= 0.5)
{
$nbrd_traits++;
}
}
elsif($x == $line-1) # bottom row
{
if($traits[$x-1][$y] >= 0.5
|| $traits[$x][$y+1] >= 0.5
|| $traits[$x][$y-1] >= 0.5
||$traits[$x-1][$y-1] >= 0.5
|| $traits[$x-1][$y+1] >= 0.5)
{
$nbrd_traits++;
}
}
elsif($y == $cols-1) # right column
{
if($traits[$x+1][$y] >= 0.5
|| $traits[$x-1][$y] >= 0.5
|| $traits[$x][$y-1] >= 0.5
||$traits[$x-1][$y-1] >= 0.5
|| $traits[$x+1][$y-1] >= 0.5)
{
$nbrd_traits++;
}
}
}
}
}
}
my $g_dispersal;
my $t_dispersal;
return ($total_genes-$nbrd_genes, $total_traits-$nbrd_traits);
}