Rotational migration
Jules Bonte
When I open population visualization however, it looks like migration seems to go from multiple subpopulations to p1 untill these subpopulations are empty.
What is weird is that even when I remove the migration block from the script, the migration continues in the same way.
I hereby provide you the whole script since I have no clue what is causing this migration.
Sincerely
Jules
// Keywords: nonWF, non-Wright-Fisher
initialize() {
initializeSLiMModelType("nonWF");
initializeSLiMOptions(keepPedigrees=T);
initializeSex("A");
defineConstant("carryingCapacities", c(400,500,600,700, 800, 900, 1000, 1100, 1200, 130));
defineConstant("frage",2);
defineConstant("lrage",7);
defineConstant("P",20);
initializeMutationType("m1", 0.5, "f", 0.0);
m1.convertToSubstitution = T;
initializeGenomicElementType("g1", m1, 1.0);
initializeGenomicElement(g1, 0, 99999);
initializeMutationRate(1e-7);
initializeRecombinationRate(1e-8);
}
reproduction(NULL, "F"){
sp = sim.subpopulations;
partnertag = individual.getValue("partner");
mate=NULL;
males=sp.individuals[sp.individuals.sex == "M" & sp.individuals.age > 3];
//if the female is single, select a single male as partner
if (partnertag == 0) {
if (males.size() > 0){
mates=males[males.getValue("partner")==0];
if (mates.size() > 0){
mate=sample(mates, 1);
individual.setValue("partner", mate.tag);
mate.setValue("partner", individual.tag);
}
}
}
else{
//if the female allready has a partner, this partner is kept
mate = subpop.subsetIndividuals(tag=partnertag);
}
random_value = runif(1);
// determine nr of offspring
if (random_value < 0.31) {
litter_size = 1;
} else if (random_value < 0.80) { // 0.31 + 0.49 = 0.80
litter_size = 2;
} else if (random_value < 0.85) { // 0.80 + 0.05 = 0.85
litter_size = 3;
} else {
litter_size = 4;
}
if (!isNULL(mate)){
if (mate.size() != 0){
offspring=subpop.addCrossed(individual, mate, count=litter_size);
for(i in seqLen(offspring.size())){
offspring[i].tag = sim.tag;
sim.tag = sim.tag + 1;
}
offspring.setValue("partner", 0);
}
}
// Disable callback
self.active = 1;
}
1 early() {
sim.addSubpop("p1", P);
p1.individuals.age = rdunif(P, min=0, max=7);
p1.individuals.tag = 1:P;
sim.tag = P+1;
for (ind in p1.individuals){
ind.setValue("partner", 0);
}
}
early() {
// life table based individual mortality
inds = sim.subpopulations.individuals;
ages = inds.age;
inds01 = inds[ages == 0];
inds12 = inds[ages == 1];
inds23 = inds[ages == 2];
inds34 = inds[ages == 3];
inds4p = inds[ages >= 4 & ages < 25];
inds25 = inds[ages >= 25];
// Mortality based on age and sex
death01_m = (runif(inds01.size()) < 0.384) & (inds01.sex == "M");
death12_m = (runif(inds12.size()) < 0.08) & (inds12.sex == "M");
death23_m = (runif(inds23.size()) < 0.048) & (inds23.sex == "M");
death34_m = (runif(inds34.size()) < 0.029) & (inds34.sex == "M");
death4p_m = (runif(inds4p.size()) < 0.15) & (inds4p.sex == "M");
death01_f = (runif(inds01.size()) < 0.284) & (inds01.sex == "F");
death12_f = (runif(inds12.size()) < 0.089) & (inds12.sex == "F");
death23_f = (runif(inds23.size()) < 0.20) & (inds23.sex == "F");
death34_f = (runif(inds34.size()) < 0.112) & (inds34.sex == "F");
death4p_f = (runif(inds4p.size()) < 0.112) & (inds4p.sex == "F");
// Kill individuals
dead_inds = (c(inds01[death01_m], inds12[death12_m],
inds23[death23_m], inds34[death34_m], inds4p[death4p_m],
inds01[death01_f], inds12[death12_f],
inds23[death23_f], inds34[death34_f], inds4p[death4p_f],
inds25));
for (ind in dead_inds) {
partnertag = ind.getValue("partner");
//catn(ind.tag +"+"+partnertag);
mate = sim.subpopulations.subsetIndividuals(tag=partnertag);
mate.setValue("partner", 0);
}
//line5 = paste(community.tick, dead_inds);
//writeFile("dead.txt", line5, append=T);
sim.killIndividuals(dead_inds);
}
15 early() {
for (spop in 2:8)
{sim.addSubpop(spop,0);
}
}
16 early (){
subpops = sim.subpopulations;
catn(subpops);
for (i in 1:7) {
migrants = p1.sampleIndividuals(15);
subpops[i].takeMigrants(migrants);
}
}
17:100 early() {
//fitness Scaling
for (i in 0:(sim.subpopulations.size() - 1)) {
subpop = sim.subpopulations[i];
Ksub = carryingCapacities[i];
if(subpop.individuals.size() > 0){
subpop.fitnessScaling = Ksub / (subpop.individualCount);
}
}
//Set Rotational migration
sourcePop = sim.subpopulations[i];
destPop = sim.subpopulations[asInteger((i + 1) % sim.subpopulations.size())];
// Select all single men
migrants = sourcePop.individuals[(sourcePop.individuals.sex == "M") & (sourcePop.individuals.getValue("partner") == 0)];
if (migrants.size() > 0) {
destPop.takeMigrants(migrants);
}
}
Ben Haller
Some observations as I look at your script:
- you're setting self.active to 1 at the end of your reproduction() callback, with the intent (as stated in the comment) of disabling the callback because you're doing "big bang" reproduction. Setting it to 1 doesn't disable it; you should set it to 0. So you're doing a lot more reproducing than you intend to. (That is masked by the fact that density-dependent fitness brings the population sizes back down before you ever see the big population size in SLiMgui, I think.)
- you're using the pseudo-parameter `subpop` in your reproduction() callback code even though you're doing "big bang" reproduction. With that reproduction strategy, you should not use `subpop` or `individual` at all, since they refer to the particular focal individual that SLiM called your code to reproduce, but you are instead reproducing *everyone*. So excise those pseudo-parameters from your code. You need to instead decide which subpop each new offspring ought to go into; probably based on the subpopulation that the parents belong to.
- your life table code could be a bit more efficient; you're drawing random numbers for individuals that aren't in the sex category you're evaluating anyway, which is a waste of computational effort. Instead of:
death01_m = (runif(inds01.size()) < 0.384) & (inds01.sex == "M");
death01_f = (runif(inds01.size()) < 0.284) & (inds01.sex == "F");
I'd suggest:
death01 = (runif(inds01.size()) < ifelse(inds01.sex == "M", 0.384, 0.284));
perhaps? Not a big deal, probably, since your population sizes are small; just a thing I noticed, and wondered how I would fix. :->
- "for (i in 0:(sim.subpopulations.size() - 1))" should be "for (seqAlong(sim.subpopulations))" – just a cleaner syntax for the same thing
- I think your rotational migration code needs a redesign. For one thing, it looks like you intend it to be inside the for loop over the subpopulations, but it is outside that loop, so it only runs once and moves individuals only from the last subpop to the first. The other problem is that it is designed to process one subpopulation at a time, so there is the potantial for one individual to migrate from p1 to p2, then p2 to p3, then p3 to p4, etc., all in a single tick, because it keeps getting chosen to migrate over and over. See recipe 15.16 for discussion of this problem. That section should also point you towards possible solutions.
Those are the issues I see with a quick skim. I hope this helps; happy modeling!
Cheers,
-B.
Benjamin C. Haller
Messer Lab
Cornell University
Jules Bonte wrote on 3/24/25 3:59 AM:
--
SLiM forward genetic simulation: http://messerlab.org/slim/
---
You received this message because you are subscribed to the Google Groups "slim-discuss" group.
To unsubscribe from this group and stop receiving emails from it, send an email to slim-discuss...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/slim-discuss/ffc7b754-a7a1-412c-a090-595be4511d2fn%40googlegroups.com.