Identical (monozygotic) twins
ryan....@gmail.com
Hello wonderful SLiM community,
Is there a way to simulate monozygotic (MZ) twins in SLiM? I could not find any recipes or examples in the manual.
Based on my experience with SLiM, I assume this would only be possible in a nonWF model. Perhaps it could be accomplished with a modifyChild() callback, or see an idea below.
In my simulation, I have something like a human reproductive pattern, with many full sibs and some half sibs. I would like to include MZ twins at some small rate, say 1%.
So for example if I just made a new offspring within a reproduction() callback:
```
child = subpop.addCrossed(mother, bio_father);
```
Could I then make a MZ twin of it like so?
```
mz_twin = subpop.addCrossed(mother, bio_father);
mz_twin.haplosomes = child.haplosomes
```
It doesn’t seem fully accommodated by SLiM, e.g. the relatedness() as calculated by SLiM I would still expect to be be 0.5 (i think). And I think I would need to log the specific twin pairs as needed. But otherwise this seems like it might work.
Is this a reasonable way to make MZ twins? Is there anything to watch out for?
Thanks
Ryan
Ben Haller
The approach you suggest will not work; you will find, I think, that trying to assign into mz_twin.haplosomes gives an error because that property is read-only (or it certainly should be! :->).
Another approach that does not work, but seems like it might, is calling subpop.addCloned(child) to make a clone of child immediately after generating it; SLiM will give an error that the parent for addCloned() must not be a newly generated offspring. (That restriction could conceivably be relaxed when tree-sequence recording is off; when tree-sequence recording is on, it gets more complicated since – I think – parents must be strictly older than their offspring according to tskit. Since SLiM tries to provide the same experience with and without tree-sequence recording, the restriction is always in force.)
A third approach that occurred to me, but that also does not work, is to set a new random number seed just before generating twin 1, generate twin 1, set the same seed again, and generate twin 2. This works in the sense that you get two identical twins; the problem is that the new mutations they acquire during meiosis will be identical mutations, BUT will be represented by distinct mutation objects, not by the same mutation objects. This would lead to incorrect fitness calculations downstream, so don't do this.
So, the correct approach is to use the addRecombinant() method (or, if you're simulating multiple chromosomes, the addMultiRecombinant() method). These methods let you make two offspring from the same parental haplosomes, with the same initial copy strand and the same crossover breakpoints; use randomizeStrands=F and choose a random initial copy strand yourself instead so that that is the same for both twins, and use the drawBreakpoints() method of Chromosome to draw breakpoints in the standard SLiM way. The recipes in sections 16.4, 16.5, and 16.9 might be of interest for understanding the details of this; there's some subtlety here since this is quite a low-level method. A tricky thing here is that even if you use one of these methods to produce two offspring that are generated by the identical meiotic process (same parental haplosomes, same initial copy strand, same crossover breakpoints) they will have independently generated mutations from their meiosis, so they will not be genetically quite identical. Apart from those new mutations they will be 100% identical by descent, but they will differ in their new mutations. There are two options for that:
- Shrug and say that that's fine; maybe you regard those mutational differences as being, as an approximation, mitotic mutations that occurred after the zygote split;
- Get a bit deeper into the weeds. You could allow the first of the two offspring that you generate to have normal SLiM-generated mutations. Then, before the second call to add[Multi]Recombinant(), set the mutation rate to zero, generate the second twin with add[Multi]Recombinant() with the same parameters, and set the mutation rate back the way it was. And then, finally, copy all of the new mutations in twin 1 (those mutations with an origin tick equal to the current tick) into twin 2 using addMutations(). Then your two twins are 100% identical, and the new mutations they share use the same mutation objects and so fitness calculations in downstream generations will be correct.
Here's a complete worked example that copies the mutations from twin1 to twin2:
initialize() {
defineConstant("MU", 1e-7); // mutation rate
defineConstant("P_MZTWIN", 0.02); // probability of an MZ twin
defineConstant("K", 500); // carrying capacity
initializeSLiMModelType("nonWF");
initializeSex();
initializeMutationType("m1", 0.5, "f", 0.0).convertToSubstitution = T;
initializeGenomicElementType("g1", m1, 1.0);
initializeGenomicElement(g1, 0, 1e7 - 1);
initializeMutationRate(MU);
initializeRecombinationRate(1e-8);
}
reproduction(NULL, "F") {
mate = subpop.sampleIndividuals(1, sex="M");
if (rdunif(1) < (1 - P_MZTWIN))
{
// make two non-twinned offspring most of the time
child1 = subpop.addCrossed(individual, mate);
child2 = subpop.addCrossed(individual, mate);
}
else
{
// make monozygotic twins with probability P_MZTWIN
ind_haplos = sample(individual.haplosomes, 2, replace=F);
mate_haplos = sample(mate.haplosomes, 2, replace=F);
breaks1 = sim.chromosome.drawBreakpoints(individual);
breaks2 = sim.chromosome.drawBreakpoints(mate);
sex = sample(c("M", "F"), 1);
tick = community.tick;
// twin1 gets new mutations as usual
twin1 = subpop.addRecombinant(ind_haplos[0], ind_haplos[1], breaks1,
mate_haplos[0], mate_haplos[1], breaks2, sex,
individual, mate, randomizeStrands=F);
// twin2 gets no new mutations
sim.chromosome.setMutationRate(0.0);
twin2 = subpop.addRecombinant(ind_haplos[0], ind_haplos[1], breaks1,
mate_haplos[0], mate_haplos[1], breaks2, sex,
individual, mate, randomizeStrands=F);
sim.chromosome.setMutationRate(MU);
// twin2 copies new mutations from twin1
twin1_muts_hap0 = twin1.haplosomes[0].mutations;
twin1_newmuts_hap0 = twin1_muts_hap0[twin1_muts_hap0.originTick == tick];
twin2.haplosomes[0].addMutations(twin1_newmuts_hap0);
twin1_muts_hap1 = twin1.haplosomes[1].mutations;
twin1_newmuts_hap1 = twin1_muts_hap1[twin1_muts_hap1.originTick == tick];
twin2.haplosomes[1].addMutations(twin1_newmuts_hap1);
}
}
1 early() {
sim.addSubpop("p1", K);
}
early() {
p1.fitnessScaling = K / p1.individualCount;
}
2000 late() { sim.outputFixedMutations(); }
Probably I'll make this a new recipe in the next version of SLiM, since it's a bit tricky. I made it a separate-sexes model to show how the sex of the offspring can be chosen and then specified as identical for both twins. A multi-chromosome version of this would be a fair bit more work, and is an exercise for the reader; but it would follow the same logic, it's just that addMultiRecombinant() is quite a complicated API since it manages so many details.
You might notice that relatedness() returns 0.5 for the twins. This is expected; relatedness() calculates relatedness based upon SLiM's pedigree-tracking information, so it is pedigree relatedness, not actual genetic relatedness. Since the twins share both parents, they are siblings and thus have relatedness() of 0.5. I wouldn't regard that as a bug; it is simply a consequence of relatedness() being pedigree relatedness. (Even for non-twin siblings the actual genetic relatedness would vary from 0.5.) You could use genetic markers or tree-sequence recording to calculate the genetic relatedness of individuals.
Speaking of tree-sequence recording, I think the copied mutations should be recorded correctly in the tree sequence, but it might look a bit funny, I'm not sure; I'd have to think about that. I think it would be apparent that the mutations actually arose in twin1 and were copied over to twin2, but I'm not 100% sure how it gets recorded in this case without doing some thinking that I don't want to do right now. :->
I hope this helps! Happy modeling!
Cheers,
-B.
Benjamin C. Haller
Messer Lab
Cornell University
--
SLiM forward genetic simulation: http://messerlab.org/slim/
Before posting, please read http://benhaller.com/slim/bugs_and_questions.html
---
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/f1e9d95f-e129-4869-87ac-890435d60169n%40googlegroups.com.