Population size causing simulation to finish before N_generations have completed

15 views
Skip to first unread message

Nkrumah Grant

unread,
Aug 12, 2021, 1:02:38 PM8/12/21
to slim-discuss

Hello, 

I am a new user of SLiM and am running into an issue simulating haploid genomes within the framework. Specifically, my SLiM simulation does not run to the generation I have defined when Ne > 130 individuals. The immature termination occurs independent of the mutation rate, genome size or number of generations I define. A reproducible example is below. It should be noted that I am running this script in the SLiM GUI version 3.6. Any help addressing this issue will be greatly appreciated.

Nkrumah


initialize()

{

defineConstant("Ne", 100);

defineConstant("Mu", 1e-9);

defineConstant("genomeSize", 1e6);

defineConstant("N_generations", 160000);

initializeMutationRate(Mu);

initializeMutationType("m1", 1.0, "f", 0.2); //neutral mutation

initializeMutationType("m2", 1.0, "f", 0.25); //beneficial mutation

initializeMutationType("m3", 1, "f", 0.25); //deleterious mutation

m1.color = "yellow";

m2.color = "green";

m3.color = "red";


m1.colorSubstitution = "grey";

m2.colorSubstitution = "green";

m3.colorSubstitution = "red";


initializeGenomicElementType("g1", c(m1), c(1));

initializeGenomicElementType("g2", c(m1,m2,m3), c(0.5,0.4,0.1)); // purifying selection

initializeGenomicElementType("g3", c(m1,m2,m3), c(0.2,0.7,0.1)); // positive selection

initializeGenomicElementType("g4", c(m1,m3), c(0.5,0.5)); // non-essential genes

g1.color = "yellow";

g2.color = "blue";

g3.color = "green";

g4.color = "red";


for (index in 0:9)

initializeGenomicElement(g1, index*1000, index*1000 + 499);

for (index in 10:19)

initializeGenomicElement(g2, index*1000, index*1000 + 499);

for (index in 20:29)

initializeGenomicElement(g3, index*1000, index*1000 + 499);

for (index in 30:39)

initializeGenomicElement(g4, index*1000, index*1000 + 499);

initializeRecombinationRate(0); // In SLiM recombination is between sister chromatids

}


// Create subpopulation 


1 late ()

{

sim.chromosome.colorSubstitution = ""; //retains color of mutation events set in the initialize callback.

sim.addSubpop("p1", asInteger(Ne));

p1.setCloningRate(1.0); // Essential for bacteria.

sim.rescheduleScriptBlock(s1, start=N_generations, end=N_generations);

}


// Remove mutation on the 2nd chromosome

// at each generation so we don't keep in memory useless mutations

modifyChild()

{

childGenome2.removeMutations();

return T;

}


// At each generation:

// - Remove fixed mutations

late()

{

// Remove fixed mutation

// (SLiM does it automatically for diploid,

// but for haploid, mutations are fixed at 0.5)

muts_m1 = sim.mutationsOfType(m1);

freqs = sim.mutationFrequencies(NULL, muts_m1);

fixed_muts_m1 = muts_m1[freqs == 0.5];

if (size(fixed_muts_m1) > 0)

sim.subpopulations.genomes.removeMutations(fixed_muts_m1, T);

muts_m2 = sim.mutationsOfType(m2);

freqs = sim.mutationFrequencies(NULL, muts_m2);

fixed_muts_m2 = muts_m2[freqs == 0.5];

if (size(fixed_muts_m2) > 0)

sim.subpopulations.genomes.removeMutations(fixed_muts_m2, T);

muts_m3 = sim.mutationsOfType(m3);

freqs = sim.mutationFrequencies(NULL, muts_m3);

fixed_muts_m3 = muts_m3[freqs == 0.5];

if (size(fixed_muts_m3) > 0)

sim.subpopulations.genomes.removeMutations(fixed_muts_m3, T);

}


// End simulation

s1 160000 late()

{

sim.outputFull("TestFullOutput.txt");

sim.simulationFinished();

}

Ben Haller

unread,
Aug 12, 2021, 5:19:52 PM8/12/21
to slim-discuss
Hi Nkrumah!

Nice catch.  This was actually caused by three different bugs, each of which was pretty harmless on its own, but when all three happened at the same time it prevented SLiM from working (and also made it fail to report the error, which was one of the three bugs :->).  I've just fixed all three.  If you use the current GitHub head version of SLiM, you should no longer experience this problem; note that chapter 2 of the SLiM manual has instructions on obtaining and building the GitHub head version, including building SLiMgui.

Since that's a pain, though, let me note that I think the bug will probably not recur if you change your modifyChild() callback slightly, to this:

modifyChild()
{
   if (size(childGenome2.mutations))
      childGenome2.removeMutations();
   return T;
}

There is nothing wrong with your original code; this tweak just avoids the bug, at least in most circumstances, by avoiding one of the three bugs that trigger it.  Alternatively, I think adding a call to "initializeSLiMOptions(mutationRuns=1);" at the top of your initialize() callback should also avoid the bug, perhaps more reliably than the first workaround, but at the cost of somewhat reduced performance for your model since it disables an internal optimization related to the bug.

Incidentally, removing mutations from the second genomes of all offspring simultaneously, with a single vectorized call, should be much faster than doing it on a per-offspring basis in modifyChild() as your model does.  Recipe 14.9 in the SLiM manual demonstrates this technique.  If there's no compelling reason for you to do it the way you're doing it, I'd recommend that you switch to the recipe's method.  This, too, would probably work around the bug, so if you do switch, and it does fix the bug for you, then you probably don't need any of the above workarounds.  :->

Good luck, and thanks for the report!  This was quite an obscure issue, and it's good to have it fixed!

Cheers,
-B.

Benjamin C. Haller
Messer Lab
Cornell University


Reply all
Reply to author
Forward
0 new messages