How can I improve the performance (in terms of computation time) of this code?
# Settings const nbloci = 100 # length of the genome const N = 100 # Number individuals in the population const nbgenerations = 100 # number of generations const mu = 1/10^5 # mutation rate const s = 0.01 # effect of a given mutation on fitness # Packages using Distributions.Binomial using Distributions.wsample # Type and Functions type Ind Genome end function makepopulation(N) Pop = Array(Ind, N) for i = 1:N Pop[i] = Ind(ones(Float64, nbloci)) end return(Pop) end function calculatefitnesses(Pop) fitnesses = Array(Float64, N) for (i,I) = enumerate(Pop) fit = 1 for l = I.Genome fit *= l end fitnesses[i] = fit end return(fitnesses) end function mutate(I) nbmut = rand(Binomial(nbloci, mu), 1) mutposs = rand(1:nbloci, nbmut) for mutpos = mutposs I.Genome[mutpos] == 1.0 ? I.Genome[mutpos] = mutfit : I.Genome[mutpos] = 1.0 end end function reproduce(Pop) fitnesses = calculatefitnesses(Pop) newPop = Array(Ind, N) for i = 1:N newPop[i] = deepcopy(wsample(Pop, fitnesses, 1)) mutate(newPop[i]) end return(newPop) end ##### MAIN ##### # Just a useful constant const mutfit = 1.0-s # Simulation Pop = makepopulation(N) # Create population for generation = 1:nbgenerations # Iterate over all generations Pop = reproduce(Pop) end
The above code simulates a population of individuals (
Ind) that all have a genome long of
nbloci (independent) loci (a locus (=sing. of loci) is a position on a chromosome).
makepopulation creates a population of clones at the beginning of the simulation. Then, iterating over all generations (
for loop) the function
reproduce is called on the population. The function
calcfitnesses who return an array of fitnesses.
reproduce sample the individuals according to
fitnesses and then mutate the sampled individuals by calling the function
mutate. After reproducing, all individuals die (no overlap of generations).
mutate, each locus has a probability
mu of mutating. In
calcfitness an individual that has for example 8 mutant alleles (allele=variant of a gene) has a fitness of
reproduce, individuals are chosen to reproduce with some probability given by their fitness.
N designates the population size and
nbgenerations obviously designates the number of generations to simulate.
My main doubts and questions
calcfitnessfunction might eventually be improved by first counting the number of mutants and then using a power function to calculate the
Ind. In such case, I may eventually want to code
True/Falseinstead of the
1.0/1.0-sin their genome. I could as well have an Array where at position
iis the fitness of on individual that has
imutations. What do you think would be the best solution?
I feel like I used lots of functions and an eventually unnecessary type. It has the advantage to make the code easier to read, but does it significantly slow down the simulations?
I am not 100% sure about the use of the
constprefix in my settings either. Would I be better off passing all these settings to each function instead of defining them as constant?
There might have a better solution than deep-copying all individuals forming the next generation. Would it be better to sample them all first and then to deep-copy only the individuals that has been sampled more than once?