diff --git a/QtSLiM/help/SLiMHelpClasses.html b/QtSLiM/help/SLiMHelpClasses.html
index b865b888..8ac23464 100644
--- a/QtSLiM/help/SLiMHelpClasses.html
+++ b/QtSLiM/help/SLiMHelpClasses.html
@@ -1050,27 +1050,30 @@
Beginning in SLiM 4.1, the count parameter dictates how many offspring will be generated (previously, exactly one offspring was generated). Each offspring is generated independently, based upon the given parameters. The returned vector contains all generated offspring, except those that were rejected by a modifyChild() callback. If all offspring are rejected, object<Individual>(0) is returned, which is a zero-length object vector of class Individual; note that this is a change in behavior from earlier versions, which would return NULL.
Note that this method is only for use in nonWF models. See addCrossed() for further general notes on the addition of new offspring individuals.
– (object<Individual>)addMultiRecombinant(object$ pattern, [Nfs$ sex = NULL], [No<Individual>$ parent1 = NULL], [No<Individual>$ parent2 = NULL], [Nl$ randomizeStrands = NULL], [integer$ count = 1], [logical$ defer = F])
-Generates a new offspring individual based upon the inheritance pattern specified by pattern, queues it for addition to the target subpopulation, and returns it. The new offspring will not be visible as a member of the target subpopulation until the end of the offspring generation tick cycle stage. The “pattern dictionary” supplied in pattern must be of class Dictionary (or a subclass of Dictionary), and more particularly, must be a dictionary of dictionaries structured in a specific way as described below. This method is a multi-chromosome version of the addRecombinant() method. For single-chromosome models, using addRecombinant() will be simpler; and it will be easier to understand this very complex method if you understand addRecombinant() first.
-The top-level “pattern dictionary” given by pattern specifies the way in which each chromosome should be handled. It can use integer keys, in which case each key is the id of a chromosome, or string keys, in which case each key is the symbol of a chromosome. In either case, a chromosome’s inheritance pattern is specified by an “inheritance dictionary” in pattern attached to such an id or symbol key. That inheritance dictionary should itself contain up to six keys, with the standard names "strand1", "strand2", "breaks1", "strand3", "strand4", and "breaks2". Any key which is missing in an inheritance dictionary is assumed to have a value of NULL. These key-value pairs are used in precisely the same way as the parameters of the same names for addRecombinant(), to produce the offspring haplosome(s) for each specified chromosome. There is some complication regarding how these six values can be used to produce results like crossing, cloning, and selfing, involving as many as four different “parents” for each chromosome; rather than repeating all of that documentation here, please see the addRecombinant() documentation for more information. When an inheritance dictionary is supplied for a particular chromosome, this method uses the six values that dictionary contains (strand1, strand2, breaks1, strand3, strand4, breaks2) in exactly the same way as addRecombinant() does; addMultiRecombinant() simply supports multiple chromosomes. In addition to that, however, addMultiRecombinant() also allows the pattern dictionary to omit the inheritance dictionaries for particular chromosomes; the behavior of addMultiRecombinant() in that case will be described below, after discussing all other aspects of the method’s implementation.
+Generates a new offspring individual based upon the inheritance pattern specified by pattern, queues it for addition to the target subpopulation, and returns it. The new offspring will not be visible as a member of the target subpopulation until the end of the offspring generation tick cycle stage. The “pattern dictionary” supplied in pattern must be of class Dictionary (or a subclass of Dictionary), and more particularly, must be a dictionary of dictionaries structured in a specific way as described below. This method is a multi-chromosome version of the addRecombinant() method. For single-chromosome models, using addRecombinant() will be simpler; and it will be easier to understand this extremely complex method if you understand addRecombinant() first.
+The top-level “pattern dictionary” given by pattern specifies the way in which each chromosome should be handled. It can use integer keys, in which case each key is the id of a chromosome, or string keys, in which case each key is the symbol of a chromosome. In either case, a chromosome’s inheritance pattern is specified by an “inheritance dictionary” in pattern attached to such a chromosome id or symbol key. That inheritance dictionary should itself contain up to six keys, with the standard names "strand1", "strand2", "breaks1", "strand3", "strand4", and "breaks2". Any key which is missing in an inheritance dictionary is assumed to have a value of NULL, and missing keys will be referred to having a value of NULL here for simplicity. These key-value pairs are used in precisely the same way as the parameters of the same names for addRecombinant(), to produce the offspring haplosome(s) for each specified chromosome. There is some complication regarding how these six values can be used to produce results like crossing, cloning, and selfing, involving as many as four different “parents” for each chromosome; rather than repeating all of that documentation here, please see the addRecombinant() documentation for more information. When an inheritance dictionary is supplied for a particular chromosome, this method uses the six values that dictionary contains (strand1, strand2, breaks1, strand3, strand4, breaks2) in exactly the same way as addRecombinant() does; addMultiRecombinant() simply supports multiple chromosomes. In addition to that, however, addMultiRecombinant() also allows the pattern dictionary to omit the inheritance dictionaries for particular chromosomes; the behavior of addMultiRecombinant() in that special case will be described below, after discussing all other aspects of the method’s implementation.
The sex parameter optionally specifies the sex of the offspring. The default value of NULL for sex specifies “default behavior”; in a non-sexual model this is the only legal value, and produces a hermaphroditic offspring. In a sexual model, the “default behavior” of NULL is that the offspring’s sex is dictated by the haplosome structure it inherits. For example, if pattern specifies that the offspring will have two non-null haplosomes for a chromosome of type "X", the sex of the offspring must therefore be female, whereas if it will have one non-null "X" haplosome and one null "X" haplosome, it must therefore be male. SLiM will scan through pattern to determine such constraints and enforce them. If the constraints implied by pattern are not self-consistent (if the offspring would have two non-null "X" haplosomes but also a non-null "Y" haplosome, for example), an error will be raised. The constraints defined by each chromosome type, as described in initializeChromosomeType(), must always be followed, even when using addMultiRecombinant(). If sex is NULL and pattern imposes no constraints upon the sex of the offspring, the offspring will be chosen as male or female with equal probability. A value of "M" or "F" for sex specifies that the offspring should be male or female, respectively. A float value from 0.0 to 1.0 for sex provides the probability that the offspring will be male; a value of 0.0 will produce a female, a value of 1.0 will produce a male, and for intermediate values SLiM will draw the sex of the offspring randomly according to the specified probability. In these cases where sex is not NULL, SLiM will first determine the sex of the individual as just described, and will then scan through pattern to confirm that it is compatible with the sex that was determined. Again, if there is a conflict an error will be raised; you cannot specify the sex of an individual to be incompatible with the haplosomes that it inherits, and if you specify a sex with a string or float value it is up to you to ensure that that is compatible with the specifications in pattern. (If you need more flexibility, you should probably not use a sexual model at all, but simply use chromosome types "A" and "H" in a non-sexual model, track the sex of individuals yourself with a tag value such as tagL0, and manipulate haplosomes during reproduction however you wish; SLiM then imposes no constraints.)
-By default, the offspring is considered to have no parents, since there may be more than two “parents” in the general case. If desired, parent1 and/or parent2 may be passed to explicitly establish particular individuals as the parents of the offspring, for purposes of pedigree tracking as well as for other purposes described below. In this case, if only one of parent1 and parent2 is non-NULL, that individual will be set as both of the parents of the offspring, mirroring the way that pedigree tracking works for other cases such as addCloned() and addSelfed(). It is not required for parent1 or parent2 to actually be a genetic parent of the offspring at all, although typically they would be.
-The randomizeStrands parameter is used to control the recombination behavior of addMultiRecombinant(). An inheritance dictionary can specify two parental strands with crossover breakpoints between them to generate one offspring strand with recombination – for example, with the "strand1", "strand2", and "breaks1" keys, as described above, but the same is true of the "strand3", "strand4", and "breaks2" keys, and the discussion that follows applies to both cases. If randomizeStrands is F, the supplied strands are used precisely as given; for example, "strand1" will be the initial copy strand when generating the first gamete to form the offspring. This mode should be used if you want explicit control over the initial copy strand; one example would be if your script is explicitly generating all four of the products of a meiosis event. If randomizeStrands is T, then if "strand1" and "strand2" are both non-NULL, 50% of the time they will be swapped, making "strand2" the initial copy strand for the first gamete instead. This mode (randomizeStrands=T) is usually the desired behavior, to avoid an inheritance bias due to a lack of randomization in the initial copy strand, so passing T for randomizeStrands is recommended unless you specifically desire otherwise. The default value of randomizeStrands is NULL in order to force either T or F to be explicitly chosen whenever it would make a difference; if it is left as NULL, an error will be raised if generation of the specified offspring involves recombination, since then SLiM needs to know whether the value is T or F. (This unconventional approach has been adopted because the default value was F prior to SLiM 5, but T is almost always the correct behavior, as explained above. To try to prevent accidental bugs, the new policy was thus adopted to force the user to explicitly choose T or F.)
-Finally, count is the number of offspring to generate using the given pattern and parameters, and defer is used for deferral of offspring generation, as described for addRecombinant(). Other details, such as which callbacks are involved, how meanParentAge for the offspring is calculated, gene conversion, spatial positioning of the offspring, and so forth, are all as described for addRecombinant().
-Constructing a well-formed pattern dictionary with inheritance dictionaries for every chromosome can be a bit complex and require many lines of code. To ease the process, see the Chromosome methods addPatternForCross(), addPatternForClone(), addPatternForNull(), and addPatternForRecombinant(), which help you to build a pattern dictionary one inheritance dictionary at a time. However, several of these methods will probably be used infrequently, because of the final aspect of addMultiRecombinant() that we have not yet discussed.
+By default, the offspring is considered to have no parents, since there may be more than two “parents” in the general case. If specifying parentage is desired, parent1 and/or parent2 may be passed explicitly; this will establish those individuals as the parents of the offspring for purposes of pedigree tracking, and for several other purposes described below. If only one of parent1 and parent2 is non-NULL, that individual will be set as both of the parents of the offspring, mirroring the way that parentage is tracked for other cases such as addCloned() and addSelfed(). It is not required for parent1 or parent2 to actually be a genetic parent of the offspring at all, although typically they would be. To benefit from the full functionality of addMultiRecombinant() as described below, it is best to supply parent1 and parent2 when possible.
+The randomizeStrands parameter is used to control the recombination behavior of addMultiRecombinant(). An inheritance dictionary can specify two parental strands with crossover breakpoints between them to generate one offspring strand with recombination – for example, with the "strand1", "strand2", and "breaks1" keys, as described above, but the same is true of the "strand3", "strand4", and "breaks2" keys, and the discussion that follows applies to both cases. If randomizeStrands is F, the supplied strands are used as given; for example, "strand1" will be the initial copy strand when generating the first gamete to form the offspring. This mode should be used if you want explicit control over the initial copy strand; one example would be if your script is explicitly generating all four of the products of a meiosis event. If randomizeStrands is T, then if "strand1" and "strand2" are both non-NULL, 50% of the time they will be swapped, making "strand2" the initial copy strand for the first gamete instead. This mode (randomizeStrands=T) is usually the desired behavior, to avoid an inheritance bias due to a lack of randomization in the initial copy strand, so passing T for randomizeStrands is recommended unless you specifically desire otherwise. The default value of randomizeStrands is NULL in order to force either T or F to be explicitly chosen whenever it would make a difference; if it is left as NULL, an error will be raised if generation of the specified offspring involves recombination, since then SLiM needs to know whether the value is T or F. (This unconventional approach has been adopted because the default value was F prior to SLiM 5, but T is almost always the correct behavior, as explained above. To try to prevent accidental bugs, this new policy was adopted to force the user to explicitly choose T or F whenever it matters.)
+The value of the meanParentAge property of the generated offspring is calculated from the mean parent age of each of its two haplosomes (whether they turn out to be null haplosomes or not). The addRecombinant() documentation provides a simple example; for addMultiRecombinant() the logic is the same, but potentially extended to more than two offspring haplosomes.
+Callbacks can be involved in offspring generation with addMultiRecombinant(), but because there are up to four strands with up to four different parents, things are a bit complicated and different from other add...() methods; the policy described here seems like the best compromise. The target subpopulation for the addMultiRecombinant() call will be used to locate applicable mutation() and modifyChild() callbacks governing the generation of the offspring individual. On the other hand, recombination() callbacks will be found based upon the subpopulations to which parent1 and parent2 belong (for reasons discussed further in drawBreakpoints()); if a parent individual is not supplied, recombination() callbacks will not be called at all when generating the corresponding offspring haplosome.
+When breakpoints are explicitly supplied to addMultiRecombinant() with breaks1 or breaks2, gene conversion tracts are not well-supported by this method; the breaks1 and breaks2 vectors provide simple crossover breakpoints, which may be used to implement crossovers or simple gene conversion tracts, but complex gene conversion tracts with heteroduplex mismatch repair are not supported in this mode of operation since there is no way to supply the relevant information. If, on the other hand, breaks1 or breaks2 is NULL when generating a haplosome with recombination, then as described above, addRecombinant() will generate breakpoints internally for that cross, and in this case, complex gene conversion tracts with heteroduplex mismatch repair are supported, since all of the necessary information is available. Similarly, if the inheritance dictionary for a given chromosome is omitted from pattern entirely and a cross between parent1 and parent2 is inferred (as discussed below), the recombination algorithm used will support gene conversion including heteroduplex mismatch repair.
+Finally, count is the number of offspring to generate using the given pattern and parameters, and defer is used for deferral of offspring generation, as described for addRecombinant(). Any other details omitted from this documentation are all as described for addRecombinant().
+Constructing a well-formed pattern dictionary with inheritance dictionaries for every chromosome can be a bit complex and require many lines of code. To ease the process, see the Chromosome methods addPatternForCross(), addPatternForClone(), addPatternForNull(), and addPatternForRecombinant(), which help you to build a pattern dictionary one inheritance dictionary at a time. However, several of these methods will probably be used infrequently, because of the final aspect of addMultiRecombinant() that we have not yet properly discussed.
As mentioned earlier, not all chromosomes need to be specified with an inheritance dictionary in pattern; whenever SLiM’s default inheritance behavior is well-defined and is desired for a given chromosome, the inheritance dictionary for that chromosome may be omitted, and will be inferred by addMultiRecombinant() automatically. This behavior makes it easy to specify a reproduction event that is, for example, like a regular biparental cross involving many chromosomes, but that uses a different reproductive pattern just for one particular chromosome that behaves in a special way. The inferred inheritance dictionary for a given chromosome is based upon the chromosome type, the values of parent1 and parent2, and the sex of the offspring (which has, by this point, been determined in all cases). The rules for this inference are actually quite simple. If both parents are specified (that is, are both non-NULL), the inferred inheritance dictionary is the same as would be produced by a call to addPatternForCross() with those two parents, given in that order, for that chromosome, for the determined offspring sex. If only one parent is specified (non-NULL), the inferred inheritance dictionary is the same as would be produced by a call to addPatternForClone() with that one parent, for that chromosome, for the determined offspring sex. If selfing is desired for the inferred inheritance dictionary, pass the same individual for both parent1 and parent2; the behavior of addPatternForCross() in that case is essentially to self the individual, as discussed in that method. If the inferred inheritance dictionary for a given chromosome is not well-defined, as discussed in the documentation for addPatternForCross() and addPatternForClone(), or if both parent1 and parent2 are NULL, an error will be raised. In such cases, the inheritance dictionary cannot be inferred, and will need to be given explicitly in pattern.
– (object<Individual>)addRecombinant(No<Haplosome>$ strand1, No<Haplosome>$ strand2, Ni breaks1, No<Haplosome>$ strand3, No<Haplosome>$ strand4, Ni breaks2, [Nfs$ sex = NULL], [No<Individual>$ parent1 = NULL], [No<Individual>$ parent2 = NULL], [Nl$ randomizeStrands = NULL], [integer$ count = 1], [logical$ defer = F])
-Generates a new offspring individual from the given parental haplosomes with the specified crossover breakpoints, queues it for addition to the target subpopulation, and returns it. The new offspring will not be visible as a member of the target subpopulation until the end of the offspring generation tick cycle stage. The target subpopulation will be used to locate applicable mutation() and modifyChild() callbacks governing the generation of the offspring individual (unlike the other addX() methods, because there are potentially up to four parental individuals to reference); recombination() callbacks will not be called by this method. This method is an advanced feature; most models will use addCrossed(), addSelfed(), or addCloned() instead. This method may only be used in single-chromosome models; in multi-chromosome models, use addMultiRecombinant(), a more general version of this method.
-This method supports several possible configurations for strand1, strand2, and breaks1 (and the same applies for strand3, strand4, and breaks2). If strand1 and strand2 are both NULL, the corresponding haplosome in the generated offspring will be empty, as from addEmpty(), with no parental haplosomes and no added mutations; in this case, breaks1 must be NULL or zero-length. If strand1 is non-NULL but strand2 is NULL, the corresponding haplosome in the generated offspring will be a clonal copy of strand1 with mutations added, as from addCloned(); in this case, breaks1 must similarly be NULL or zero-length. If strand1 and strand2 are both non-NULL, the corresponding haplosome in the generated offspring will result from recombination between strand1 and strand2 with mutations added, as from addCrossed(), with strand1 being the initial copy strand; copying will switch between strands at each breakpoint in breaks1, which must be non-NULL but need not be sorted or uniqued (SLiM will sort and unique the supplied breakpoints internally). (It is not currently legal for strand1 to be NULL and strand2 non-NULL; that variant may be assigned some meaning in future.) Again, this discussion applies equally to strand3, strand4, and breaks2, mutatis mutandis. Note that when new mutations are generated by addRecombinant(), their subpopID property will be the id of the offspring’s subpopulation, since the parental subpopulation is ambiguous in the general case; this behavior differs from the other add...() methods.
-The sex parameter is interpreted exactly as in addCrossed(); see that method for discussion. If the offspring sex is specified in any way (i.e., if sex is non-NULL), the strands provided must be compatible with the sex chosen. If the offspring sex is not specified (i.e., if sex is NULL), the sex will be inferred from the strands provided where possible (when modeling a chromosome with a sex-specific inheritance pattern), or will be chosen randomly otherwise (e.g., when modeling autosomes); it will not be inferred from the sex of the individuals possessing the parental strands, even when the reproductive mode is essentially clonal from a single parent, since such inference would be ambiguous in the general case.
-By default, the offspring is considered to have no parents for the purposes of pedigree tracking, since there may be more than two “parents” in the general case. If pedigree tracking of parentage is desired, parent1 and/or parent2 may be passed to explicitly establish particular individuals as the parents of the offspring for purposes of pedigree tracking. In this case, if only one of parent1 and parent2 is non-NULL, that individual will be set as both of the parents of the offspring, mirroring the way that parentage is tracked for other cases such as addCloned() and addSelfed(). It is not required for parent1 or parent2 to actually be a genetic parent of the offspring at all, although typically they would be.
-The randomizeStrands parameter is used to control the recombination behavior of addRecombinant(). Two parental strands with crossover breakpoints between them can be specified as described above to generate one offspring strand with recombination – for example, with the strand1, strand2, and breaks1 parameters, as described above, but the same is true of the strand3, strand4, and breaks2 parameters, and the discussion that follows applies to both cases. If randomizeStrands is F, the supplied strands are used precisely as given; for example, strand1 will be the initial copy strand when generating the first gamete to form the offspring. This mode should be used if you want explicit control over the initial copy strand; one example would be if your script is explicitly generating all four of the products of a meiosis event. If randomizeStrands is T, then if strand1 and strand2 are both non-NULL, 50% of the time they will be swapped, making strand2 the initial copy strand for the first gamete instead. This mode (randomizeStrands=T) is usually the desired behavior, to avoid an inheritance bias due to a lack of randomization in the initial copy strand, so passing T for randomizeStrands is recommended unless you specifically desire otherwise. The default value of randomizeStrands is NULL in order to force either T or F to be explicitly chosen whenever it would make a difference; if it is left as NULL, an error will be raised if generation of the specified offspring involves recombination, since then SLiM needs to know whether the value is T or F. (This unconventional approach has been adopted because the default value was F prior to SLiM 5, but T is almost always the correct behavior, as explained above. To try to prevent accidental bugs, the new policy was thus adopted to force the user to explicitly choose T or F.)
-If randomizeStrands is F (the default), strand1 will be the initial copy strand when generating the first gamete to form the offspring, and strand3 will be the initial copy strand when generating the second gamete. If randomizeStrands is T, then if strand1 and strand2 are both non-NULL, 50% of the time they will be swapped, making strand2 the initial copy strand for the first gamete; and similarly, if strand3 and strand4 are both non-NULL, 50% of the time they will be swapped, making strand4 the initial copy strand for the second gamete. This is probably usually the desired behavior, to avoid an inheritance bias due to a lack of randomization in the initial copy strand, so passing T for randomizeStrands is recommended unless you specifically desire otherwise. It is not the default behavior only for reasons of backward compatibility.
-These semantics allow several uses for addRecombinant(). When all strands are non-NULL, it is similar to addCrossed() except that the recombination breakpoints are specified explicitly, allowing very precise offspring generation without having to override SLiM’s breakpoint generation with a recombination() callback. When only strand1 and strand3 are supplied, it is very similar to addCloned(), creating a clonal offspring, except that the two parental haplosomes need not belong to the same individual (whatever that might mean biologically). Supplying only strand1 is useful for modeling clonally reproducing haploids; the second haplosome of every offspring will be kept empty and will not receive new mutations. For a model of clonally reproducing haploids that undergo horizontal gene transfer (HGT), supplying only strand1 and strand2 will allow HGT from strand2 to replace segments of an otherwise clonal copy of strand1, while the second haplosome of the generated offspring will again be kept empty; this could be useful for modeling bacterial conjugation, for example. Other variations are also possible.
-The value of the meanParentAge property of the generated offspring is calculated from the mean parent age of each of its two haplosomes (whether they turn out to be null haplosomes or not); that may be an average of two values (if both offspring haplosomes have at least one parent), a single value (if one offspring haplosome has no parent), or no values (if both offspring haplosomes have no parent, in which case 0.0 results). The mean parent age of a given offspring haplosome is the mean of the ages of the parents of the two strands used to generate that offspring haplosome; if one strand is NULL then the mean parent age for that offspring haplosome is the age of the parent of the non-NULL strand, while if both strands are NULL then that offspring haplosome is parentless and is not used in the final calculation. In other words, if one offspring haplosome has two parents with ages A and B, and the other offspring haplosome has one parent with age C, the meanParentAge of the offspring will be (A+B+C+C) / 4, not (A+B+C) / 3.
-Note that gene conversion tracts are not explicitly supported by this method; the breaks vectors provide crossover breakpoints, which may be used to implement crossovers or simple gene conversion tracts. There is no way to specify complex gene conversion tracts with heteroduplex mismatch repair.
+Generates a new offspring individual from the given parental haplosomes with the specified crossover breakpoints, queues it for addition to the target subpopulation, and returns it. The new offspring will not be visible as a member of the target subpopulation until the end of the offspring generation tick cycle stage. This method is an advanced feature; most models will use addCrossed(), addSelfed(), or addCloned() instead. This method may only be used in single-chromosome models; in multi-chromosome models, use addMultiRecombinant(), a more general version of addRecombinant().
+This method supports several possible configurations for strand1, strand2, and breaks1 (and the same applies for strand3, strand4, and breaks2). If strand1 and strand2 are both NULL, the corresponding haplosome in the generated offspring will be a null haplosome; in this case, breaks1 must be NULL or zero-length. If strand1 is non-NULL but strand2 is NULL, the corresponding haplosome in the generated offspring will be a clonal copy of strand1 with mutations added, as from addCloned(); in this case, breaks1 must again be NULL or zero-length. If strand1 and strand2 are both non-NULL, the corresponding haplosome in the generated offspring will result from recombination between strand1 and strand2 with mutations added, as from addCrossed(), with strand1 being the initial copy strand by default (but see below). Copying will switch between strands at each crossover breakpoint. Breakpoints may be supplied in breaks1, which need not be sorted or uniqued (SLiM will sort and unique the supplied breakpoints internally). Alternatively, breaks1 may be NULL, which requests that addRecombinant() draw breakpoints automatically to recombine strand1 and strand2, following SLiM’s usual breakpoint-drawing algorithm. (If you do not want any breakpoints, pass integer(0), a zero-length integer vector, for breaks1.) Finally, it is not currently legal for strand1 to be NULL and strand2 non-NULL; that variant may be assigned some meaning in future. Again, this discussion applies equally to strand3, strand4, and breaks2, mutatis mutandis. Null haplosomes may never be passed as any of the four parental strands; pass NULL, not a null haplosome, if that strand is not inherited from. When modeling a chromosome that is intrinsically haploid, such as the Y, NULL must be passed for strand3, strand4, and breaks2; you cannot supply genetic information for an offspring haplosome that will not exist. Note that when new mutations are generated by addRecombinant(), their subpopID property will be the id of the offspring’s subpopulation, since the parental subpopulation is ambiguous in the general case; this behavior differs from the other add...() methods.
+These semantics allow several uses for addRecombinant(). When all strands are non-NULL, it is similar to addCrossed() except that the recombination breakpoints can be specified explicitly, allowing very precise offspring generation without having to override SLiM’s breakpoint generation with a recombination() callback. When only strand1 and strand3 are supplied, it is very similar to addCloned(), creating a clonal offspring, except that the two parental haplosomes need not belong to the same individual (whatever that might mean biologically). Supplying only strand1 is useful for modeling clonally reproducing haploids, or any chromosome type that is intrinsically haploid, such as the Y chromosome. For a model of clonally reproducing haploids that undergo horizontal gene transfer (HGT), supplying only strand1 and strand2 will allow HGT from strand2 to replace segments of an otherwise clonal copy of strand1, while the second haplosome of the generated offspring will be a null haplosome; this could be useful for modeling bacterial conjugation, for example. Other variations are also possible.
+The sex parameter optionally specifies the sex of the offspring. The default value of NULL for sex specifies “default behavior”; in a non-sexual model this is the only legal value, and produces a hermaphroditic offspring. In a sexual model, the “default behavior” of NULL is that the offspring’s sex is dictated by the haplosome structure it inherits. For example, if the supplied strands indicate that the offspring will have two non-null haplosomes for a chromosome of type "X", the sex of the offspring must therefore be female, whereas if it will have one non-null "X" haplosome and one null "X" haplosome, it must therefore be male. SLiM will examine the supplied strands to determine such constraints and enforce them. The constraints defined by each chromosome type, as described in initializeChromosomeType(), must always be followed, even when using addRecombinant(). If sex is NULL and the sex of the offspring is unconstrained, the offspring will be chosen as male or female with equal probability. A value of "M" or "F" for sex specifies that the offspring should be male or female, respectively. A float value from 0.0 to 1.0 for sex provides the probability that the offspring will be male; a value of 0.0 will produce a female, a value of 1.0 will produce a male, and for intermediate values SLiM will draw the sex of the offspring randomly according to the specified probability. In these cases where sex is not NULL, SLiM will first determine the sex of the individual as just described, and will then examine the supplied strands to confirm that it is compatible with the sex that was determined. Again, if there is a conflict an error will be raised; you cannot specify the sex of an individual to be incompatible with the haplosomes that it inherits, and if you specify a sex with a string or float value it is up to you to ensure that that is compatible with the supplied strands. (If you need more flexibility, you should probably not use a sexual model at all, but simply use chromosome type "A" or "H" in a non-sexual model, track the sex of individuals yourself with a tag value such as tagL0, and manipulate haplosomes during reproduction however you wish; SLiM then imposes no constraints.)
+By default, the offspring is considered to have no parents, since there may be more than two “parents” in the general case. If specifying parentage is desired, parent1 and/or parent2 may be passed to explicitly; this will establish those individuals as the parents of the offspring for purposes of pedigree tracking, and for several other purposes described below. If only one of parent1 and parent2 is non-NULL, that individual will be set as both of the parents of the offspring, mirroring the way that parentage is tracked for other cases such as addCloned() and addSelfed(). It is not required for parent1 or parent2 to actually be a genetic parent of the offspring at all, although typically they would be. To benefit from the full functionality of addRecombinant() as described below, it is best to supply parent1 and parent2 when possible.
+The randomizeStrands parameter is used to control the recombination behavior of addRecombinant(). As described above, two parental strands with crossover breakpoints between them can be specified to generate one offspring strand with recombination – for example, with the strand1, strand2, and breaks1 parameters, but the same is true of the strand3, strand4, and breaks2 parameters, and the discussion that follows applies to both cases. If randomizeStrands is F, the supplied strands are used as given; for example, strand1 will be the initial copy strand when generating the first gamete to form the offspring. This mode should be used if you want explicit control over the initial copy strand; one example would be if your script is explicitly generating all four of the products of a meiosis event. If randomizeStrands is T, then if strand1 and strand2 are both non-NULL, 50% of the time they will be swapped, making strand2 the initial copy strand for the first gamete instead. This mode (randomizeStrands=T) is usually the desired behavior, to avoid an inheritance bias due to a lack of randomization in the initial copy strand, so passing T for randomizeStrands is recommended unless you specifically desire otherwise. The default value of randomizeStrands is NULL in order to force either T or F to be explicitly chosen whenever it would make a difference; if it is left as NULL, an error will be raised if generation of the specified offspring involves recombination, since then SLiM needs to know whether the value is T or F. (This unconventional approach has been adopted because the default value was F prior to SLiM 5, but T is almost always the correct behavior, as explained above. To try to prevent accidental bugs, this new policy was adopted to force the user to explicitly choose T or F whenever it matters.)
+The value of the meanParentAge property of the generated offspring is calculated from the mean parent age of each of its two haplosomes (whether they turn out to be null haplosomes or not); that may be an average of two values (if both offspring haplosomes have at least one parent), a single value (if one offspring haplosome has no parent), or no values (if both offspring haplosomes have no parent, in which case 0.0 results). The mean parent age of a given offspring haplosome is the mean of the ages of the parents of the two strands used to generate that offspring haplosome; if one strand is NULL then the mean parent age for that offspring haplosome is the age of the parent of the non-NULL strand, while if both strands are NULL then that offspring haplosome is parentless and is not used in the final calculation. In other words, if one offspring haplosome has two parents with ages A and B, and the other offspring haplosome has one parent with age C, the meanParentAge of the offspring will be (A+B+C+C) / 4, or equivalently, ((A+B)/2 + C) / 2, not (A+B+C) / 3.
+Callbacks can be involved in offspring generation with addRecombinant(), but because there are up to four strands with up to four different parents, things are a bit complicated and different from other add...() methods; the policy described here seems like the best compromise. The target subpopulation for the addRecombinant() call will be used to locate applicable mutation() and modifyChild() callbacks governing the generation of the offspring individual. On the other hand, recombination() callbacks will be found based upon the subpopulations to which parent1 and parent2 belong (for reasons discussed further in drawBreakpoints()); if a parent individual is not supplied, recombination() callbacks will not be called at all when generating the corresponding offspring haplosome.
+When breakpoints are explicitly supplied to addRecombinant() with breaks1 or breaks2, gene conversion tracts are not well-supported by this method; the breaks1 and breaks2 vectors provide simple crossover breakpoints, which may be used to implement crossovers or simple gene conversion tracts, but complex gene conversion tracts with heteroduplex mismatch repair are not supported in this mode of operation since there is no way to supply the relevant information. If, on the other hand, breaks1 or breaks2 is NULL when generating a haplosome with recombination, then as described above, addRecombinant() will generate breakpoints internally for that cross, and in this case, complex gene conversion tracts with heteroduplex mismatch repair are supported, since all of the necessary information is available.
Beginning in SLiM 4.1, the count parameter dictates how many offspring will be generated (previously, exactly one offspring was generated). Each offspring is generated independently, based upon the given parameters. The returned vector contains all generated offspring, except those that were rejected by a modifyChild() callback. If all offspring are rejected, object<Individual>(0) is returned, which is a zero-length object vector of class Individual; note that this is a change in behavior from earlier versions, which would return NULL.
Beginning in SLiM 4.1, passing T for defer will defer the generation of the haplosomes of the produced offspring until the end of the reproduction phase. Haplosome generation can only be deferred if there are no active mutation() callbacks; otherwise, an error will result. Furthermore, when haplosome generation is deferred the mutations of the haplosomes of the generated offspring may not be accessed until reproduction is complete (whether from a modifyChild() callback or otherwise). There is little or no advantage to deferring haplosome generation when running single-threaded; in that case, the default of F for defer is generally preferable since it has fewer restrictions. When running multi-threaded, deferring haplosome generation allows that task to be done in parallel (which is the reason this option exists).
-Also beginning in SLiM 4.1, in spatial models the spatial position of the offspring will be inherited (i.e., copied) from parent1; more specifically, the x property will be inherited in all spatial models (1D/2D/3D), the y property in 2D/3D models, and the z property in 3D models. Properties not inherited will be left uninitialized, as they were prior to SLiM 4.1. The parent’s spatial position is probably not desirable in itself; the intention here is to make it easy to model the natal dispersal of all the new offspring for a given tick with a single vectorized call to deviatePositions() / pointDeviated(). If parent1 is NULL (the default), parent2 will be used; if it is also NULL, no spatial position will be inherited.
+Also beginning in SLiM 4.1, in spatial models the spatial position of the offspring will be inherited (i.e., copied) from parent1; more specifically, the x property will be inherited in all spatial models (1D/2D/3D), the y property in 2D/3D models, and the z property in 3D models. Properties not inherited will be left uninitialized, as they were prior to SLiM 4.1. The parent’s spatial position is probably not desirable in itself; the intention here is to make it easy to model the natal dispersal of all the new offspring for a given tick with a single vectorized call to deviatePositions() / pointDeviated(). If parent1 is NULL, parent2 will be used; if it is also NULL, no spatial position will be inherited.
Note that this method is only for use in nonWF models. See addCrossed() for further general notes on the addition of new offspring individuals.
– (object<Individual>)addSelfed(object<Individual>$ parent, [integer$ count = 1], [logical$ defer = F])
Generates a new offspring individual from the given parent by selfing, queues it for addition to the target subpopulation, and returns it. The new offspring will not be visible as a member of the target subpopulation until the end of the offspring generation tick cycle stage. The subpopulation of parent will be used to locate applicable mutation(), recombination(), and modifyChild() callbacks governing the generation of the offspring individual.
diff --git a/SLiMgui/SLiMHelpClasses.rtf b/SLiMgui/SLiMHelpClasses.rtf
index 4476ff0f..b3e3d613 100644
--- a/SLiMgui/SLiMHelpClasses.rtf
+++ b/SLiMgui/SLiMHelpClasses.rtf
@@ -9547,7 +9547,7 @@ Note that this method is only for use in nonWF models. See
\f3\fs18 addRecombinant()
\f4\fs20 method. For single-chromosome models, using
\f3\fs18 addRecombinant()
-\f4\fs20 will be simpler; and it will be easier to understand this very complex method if you understand
+\f4\fs20 will be simpler; and it will be easier to understand this extremely complex method if you understand
\f3\fs18 addRecombinant()
\f4\fs20 first.\
The top-level \'93pattern dictionary\'94 given by
@@ -9562,7 +9562,7 @@ The top-level \'93pattern dictionary\'94 given by
\f3\fs18 symbol
\f4\fs20 of a chromosome. In either case, a chromosome\'92s inheritance pattern is specified by an \'93inheritance dictionary\'94 in
\f3\fs18 pattern
-\f4\fs20 attached to such an
+\f4\fs20 attached to such a chromosome
\f3\fs18 id
\f4\fs20 or
\f3\fs18 symbol
@@ -9580,7 +9580,9 @@ The top-level \'93pattern dictionary\'94 given by
\f3\fs18 "breaks2"
\f4\fs20 . Any key which is missing in an inheritance dictionary is assumed to have a value of
\f3\fs18 NULL
-\f4\fs20 . These key-value pairs are used in precisely the same way as the parameters of the same names for
+\f4\fs20 , and missing keys will be referred to having a value of
+\f3\fs18 NULL
+\f4\fs20 here for simplicity. These key-value pairs are used in precisely the same way as the parameters of the same names for
\f3\fs18 addRecombinant()
\f4\fs20 , to produce the offspring haplosome(s) for each specified chromosome. There is some complication regarding how these six values can be used to produce results like crossing, cloning, and selfing, involving as many as four different \'93parents\'94 for each chromosome; rather than repeating all of that documentation here, please see the
\f3\fs18 addRecombinant()
@@ -9604,7 +9606,7 @@ The top-level \'93pattern dictionary\'94 given by
\f3\fs18 addMultiRecombinant()
\f4\fs20 also allows the pattern dictionary to omit the inheritance dictionaries for particular chromosomes; the behavior of
\f3\fs18 addMultiRecombinant()
-\f4\fs20 in that case will be described below, after discussing all other aspects of the method\'92s implementation.\
+\f4\fs20 in that special case will be described below, after discussing all other aspects of the method\'92s implementation.\
The
\f3\fs18 sex
\f4\fs20 parameter optionally specifies the sex of the offspring. The default value of
@@ -9676,11 +9678,11 @@ The
\f4\fs20 in a non-sexual model, track the sex of individuals yourself with a tag value such as
\f3\fs18 tagL0
\f4\fs20 , and manipulate haplosomes during reproduction however you wish; SLiM then imposes no constraints.)\
-By default, the offspring is considered to have no parents, since there may be more than two \'93parents\'94 in the general case. If desired,
+By default, the offspring is considered to have no parents, since there may be more than two \'93parents\'94 in the general case. If specifying parentage is desired,
\f3\fs18 parent1
\f4\fs20 and/or
\f3\fs18 parent2
-\f4\fs20 may be passed to explicitly establish particular individuals as the parents of the offspring, for purposes of pedigree tracking as well as for other purposes described below. In this case, if only one of
+\f4\fs20 may be passed explicitly; this will establish those individuals as the parents of the offspring for purposes of pedigree tracking, and for several other purposes described below. If only one of
\f3\fs18 parent1
\f4\fs20 and
\f3\fs18 parent2
@@ -9688,7 +9690,7 @@ By default, the offspring is considered to have no parents, since there may be m
\f3\fs18 NULL
\f4\fs20 , that individual will be set as
\f1\i both
-\f4\i0 of the parents of the offspring, mirroring the way that pedigree tracking works for other cases such as
+\f4\i0 of the parents of the offspring, mirroring the way that parentage is tracked for other cases such as
\f3\fs18 addCloned()
\f4\fs20 and
\f3\fs18 addSelfed()
@@ -9696,7 +9698,13 @@ By default, the offspring is considered to have no parents, since there may be m
\f3\fs18 parent1
\f4\fs20 or
\f3\fs18 parent2
-\f4\fs20 to actually be a genetic parent of the offspring at all, although typically they would be.\
+\f4\fs20 to actually be a genetic parent of the offspring at all, although typically they would be. To benefit from the full functionality of
+\f3\fs18 addMultiRecombinant()
+\f4\fs20 as described below, it is best to supply
+\f3\fs18 parent1
+\f4\fs20 and
+\f3\fs18 parent2
+\f4\fs20 when possible.\
The
\f3\fs18 randomizeStrands
\f4\fs20 parameter is used to control the recombination behavior of
@@ -9717,7 +9725,7 @@ The
\f3\fs18 randomizeStrands
\f4\fs20 is
\f3\fs18 F
-\f4\fs20 , the supplied strands are used precisely as given; for example,
+\f4\fs20 , the supplied strands are used as given; for example,
\f3\fs18 "strand1"
\f4\fs20 will be the initial copy strand when generating the first gamete to form the offspring. This mode should be used if you want explicit control over the initial copy strand; one example would be if your script is explicitly generating all four of the products of a meiosis event. If
\f3\fs18 randomizeStrands
@@ -9755,20 +9763,71 @@ The
\f3\fs18 F
\f4\fs20 prior to SLiM 5, but
\f3\fs18 T
-\f4\fs20 is almost always the correct behavior, as explained above. To try to prevent accidental bugs, the new policy was thus adopted to force the user to explicitly choose
+\f4\fs20 is almost always the correct behavior, as explained above. To try to prevent accidental bugs, this new policy was adopted to force the user to explicitly choose
\f3\fs18 T
\f4\fs20 or
\f3\fs18 F
-\f4\fs20 .)\
+\f4\fs20 whenever it matters.)\
+The value of the
+\f3\fs18 meanParentAge
+\f4\fs20 property of the generated offspring is calculated from the mean parent age of each of its two haplosomes (whether they turn out to be null haplosomes or not). The
+\f3\fs18 addRecombinant()
+\f4\fs20 documentation provides a simple example; for
+\f3\fs18 addMultiRecombinant()
+\f4\fs20 the logic is the same, but potentially extended to more than two offspring haplosomes.\
+Callbacks can be involved in offspring generation with
+\f3\fs18 addMultiRecombinant()
+\f4\fs20 , but because there are up to four strands with up to four different parents, things are a bit complicated and different from other
+\f3\fs18 add...()
+\f4\fs20 methods; the policy described here seems like the best compromise. The target subpopulation for the
+\f3\fs18 addMultiRecombinant()
+\f4\fs20 call will be used to locate applicable
+\f3\fs18 mutation()
+\f4\fs20 and
+\f3\fs18 modifyChild()
+\f4\fs20 callbacks governing the generation of the offspring individual. On the other hand,
+\f3\fs18 recombination()
+\f4\fs20 callbacks will be found based upon the subpopulations to which
+\f3\fs18 parent1
+\f4\fs20 and
+\f3\fs18 parent2
+\f4\fs20 belong (for reasons discussed further in
+\f3\fs18 drawBreakpoints()
+\f4\fs20 ); if a parent individual is not supplied,
+\f3\fs18 recombination()
+\f4\fs20 callbacks will not be called at all when generating the corresponding offspring haplosome.\
+When breakpoints are explicitly supplied to
+\f3\fs18 addMultiRecombinant()
+\f4\fs20 with
+\f3\fs18 breaks1
+\f4\fs20 or
+\f3\fs18 breaks2
+\f4\fs20 , gene conversion tracts are not well-supported by this method; the
+\f3\fs18 breaks1
+\f4\fs20 and
+\f3\fs18 breaks2
+\f4\fs20 vectors provide simple crossover breakpoints, which may be used to implement crossovers or simple gene conversion tracts, but complex gene conversion tracts with heteroduplex mismatch repair are not supported in this mode of operation since there is no way to supply the relevant information. If, on the other hand,
+\f3\fs18 breaks1
+\f4\fs20 or
+\f3\fs18 breaks2
+\f4\fs20 is
+\f3\fs18 NULL
+\f4\fs20 when generating a haplosome with recombination, then as described above,
+\f3\fs18 addRecombinant()
+\f4\fs20 will generate breakpoints internally for that cross, and in this case, complex gene conversion tracts with heteroduplex mismatch repair are supported, since all of the necessary information is available. Similarly, if the inheritance dictionary for a given chromosome is omitted from
+\f3\fs18 pattern
+\f4\fs20 entirely and a cross between
+\f3\fs18 parent1
+\f4\fs20 and
+\f3\fs18 parent2
+\f4\fs20 is inferred (as discussed below), the recombination algorithm used will support gene conversion including heteroduplex mismatch repair.\
Finally,
\f3\fs18 count
\f4\fs20 is the number of offspring to generate using the given pattern and parameters, and
\f3\fs18 defer
\f4\fs20 is used for deferral of offspring generation, as described for
\f3\fs18 addRecombinant()
-\f4\fs20 . Other details, such as which callbacks are involved, how
-\f3\fs18 meanParentAge
-\f4\fs20 for the offspring is calculated, gene conversion, spatial positioning of the offspring, and so forth, are all as described for
+\f4\fs20 . Any other details omitted from this documentation are all as described for
\f3\fs18 addRecombinant()
\f4\fs20 .\
Constructing a well-formed pattern dictionary with inheritance dictionaries for every chromosome can be a bit complex and require many lines of code. To ease the process, see the
@@ -9783,7 +9842,7 @@ Constructing a well-formed pattern dictionary with inheritance dictionaries for
\f3\fs18 addPatternForRecombinant()
\f4\fs20 , which help you to build a pattern dictionary one inheritance dictionary at a time. However, several of these methods will probably be used infrequently, because of the final aspect of
\f3\fs18 addMultiRecombinant()
-\f4\fs20 that we have not yet discussed.\
+\f4\fs20 that we have not yet properly discussed.\
As mentioned earlier, not all chromosomes need to be specified with an inheritance dictionary in
\f3\fs18 pattern
\f4\fs20 ; whenever SLiM\'92s default inheritance behavior is well-defined and is desired for a given chromosome, the inheritance dictionary for that chromosome may be omitted, and will be inferred by
@@ -9824,15 +9883,7 @@ As mentioned earlier, not all chromosomes need to be specified with an inheritan
\f3\fs18 \cf2 \'96\'a0(object)addRecombinant(No$\'a0strand1, No$\'a0strand2, Ni\'a0breaks1, No$\'a0strand3, No$\'a0strand4, Ni\'a0breaks2, [Nfs$\'a0sex\'a0=\'a0NULL], [No$\'a0parent1\'a0=\'a0NULL], [No$\'a0parent2\'a0=\'a0NULL], [Nl$\'a0randomizeStrands\'a0=\'a0NULL], [integer$\'a0count\'a0=\'a01], [logical$\'a0defer\'a0=\'a0F])\
\pard\pardeftab720\li547\ri720\sb60\sa60\partightenfactor0
-\f4\fs20 \cf2 Generates a new offspring individual from the given parental haplosomes with the specified crossover breakpoints, queues it for addition to the target subpopulation, and returns it. The new offspring will not be visible as a member of the target subpopulation until the end of the offspring generation tick cycle stage. The target subpopulation will be used to locate applicable
-\f3\fs18 mutation()
-\f4\fs20 and
-\f3\fs18 modifyChild()
-\f4\fs20 callbacks governing the generation of the offspring individual (unlike the other
-\f3\fs18 addX()
-\f4\fs20 methods, because there are potentially up to four parental individuals to reference);
-\f3\fs18 recombination()
-\f4\fs20 callbacks will not be called by this method. This method is an advanced feature; most models will use
+\f4\fs20 \cf2 Generates a new offspring individual from the given parental haplosomes with the specified crossover breakpoints, queues it for addition to the target subpopulation, and returns it. The new offspring will not be visible as a member of the target subpopulation until the end of the offspring generation tick cycle stage. This method is an advanced feature; most models will use
\f3\fs18 addCrossed()
\f4\fs20 ,
\f3\fs18 addSelfed()
@@ -9840,7 +9891,9 @@ As mentioned earlier, not all chromosomes need to be specified with an inheritan
\f3\fs18 addCloned()
\f4\fs20 instead. This method may only be used in single-chromosome models; in multi-chromosome models, use
\f3\fs18 addMultiRecombinant()
-\f4\fs20 , a more general version of this method.\
+\f4\fs20 , a more general version of
+\f3\fs18 addRecombinant()
+\f4\fs20 .\
This method supports several possible configurations for
\f3\fs18 strand1
\f4\fs20 ,
@@ -9859,9 +9912,7 @@ This method supports several possible configurations for
\f3\fs18 strand2
\f4\fs20 are both
\f3\fs18 NULL
-\f4\fs20 , the corresponding haplosome in the generated offspring will be empty, as from
-\f3\fs18 addEmpty()
-\f4\fs20 , with no parental haplosomes and no added mutations; in this case,
+\f4\fs20 , the corresponding haplosome in the generated offspring will be a null haplosome; in this case,
\f3\fs18 breaks1
\f4\fs20 must be
\f3\fs18 NULL
@@ -9879,7 +9930,7 @@ This method supports several possible configurations for
\f3\fs18 addCloned()
\f4\fs20 ; in this case,
\f3\fs18 breaks1
-\f4\fs20 must similarly be
+\f4\fs20 must again be
\f3\fs18 NULL
\f4\fs20 or zero-length. If
\f3\fs18 strand1
@@ -9895,11 +9946,25 @@ This method supports several possible configurations for
\f3\fs18 addCrossed()
\f4\fs20 , with
\f3\fs18 strand1
-\f4\fs20 being the initial copy strand; copying will switch between strands at each breakpoint in
+\f4\fs20 being the initial copy strand by default (but see below). Copying will switch between strands at each crossover breakpoint. Breakpoints may be supplied in
+\f3\fs18 breaks1
+\f4\fs20 , which need not be sorted or uniqued (SLiM will sort and unique the supplied breakpoints internally). Alternatively,
\f3\fs18 breaks1
-\f4\fs20 , which must be non-
+\f4\fs20 may be
\f3\fs18 NULL
-\f4\fs20 but need not be sorted or uniqued (SLiM will sort and unique the supplied breakpoints internally). (It is not currently legal for
+\f4\fs20 , which requests that
+\f3\fs18 addRecombinant()
+\f4\fs20 draw breakpoints automatically to recombine
+\f3\fs18 strand1
+\f4\fs20 and
+\f3\fs18 strand2
+\f4\fs20 , following SLiM\'92s usual breakpoint-drawing algorithm. (If you do not want any breakpoints, pass
+\f3\fs18 integer(0)
+\f4\fs20 , a zero-length
+\f3\fs18 integer
+\f4\fs20 vector, for
+\f3\fs18 breaks1
+\f4\fs20 .) Finally, it is not currently legal for
\f3\fs18 strand1
\f4\fs20 to be
\f3\fs18 NULL
@@ -9907,7 +9972,7 @@ This method supports several possible configurations for
\f3\fs18 strand2
\f4\fs20 non-
\f3\fs18 NULL
-\f4\fs20 ; that variant may be assigned some meaning in future.) Again, this discussion applies equally to
+\f4\fs20 ; that variant may be assigned some meaning in future. Again, this discussion applies equally to
\f3\fs18 strand3
\f4\fs20 ,
\f3\fs18 strand4
@@ -9915,7 +9980,17 @@ This method supports several possible configurations for
\f3\fs18 breaks2
\f4\fs20 ,
\f1\i mutatis mutandis
-\f4\i0 . Note that when new mutations are generated by
+\f4\i0 . Null haplosomes may never be passed as any of the four parental strands; pass
+\f3\fs18 NULL
+\f4\fs20 , not a null haplosome, if that strand is not inherited from. When modeling a chromosome that is intrinsically haploid, such as the Y,
+\f3\fs18 NULL
+\f4\fs20 must be passed for
+\f3\fs18 strand3
+\f4\fs20 ,
+\f3\fs18 strand4
+\f4\fs20 , and
+\f3\fs18 breaks2
+\f4\fs20 ; you cannot supply genetic information for an offspring haplosome that will not exist. Note that when new mutations are generated by
\f3\fs18 addRecombinant()
\f4\fs20 , their
\f3\fs18 subpopID
@@ -9924,26 +9999,91 @@ This method supports several possible configurations for
\f4\fs20 of the offspring\'92s subpopulation, since the parental subpopulation is ambiguous in the general case; this behavior differs from the other
\f3\fs18 add...()
\f4\fs20 methods.\
+These semantics allow several uses for
+\f3\fs18 addRecombinant()
+\f4\fs20 . When all strands are non-
+\f3\fs18 NULL
+\f4\fs20 , it is similar to
+\f3\fs18 addCrossed()
+\f4\fs20 except that the recombination breakpoints can be specified explicitly, allowing very precise offspring generation without having to override SLiM\'92s breakpoint generation with a
+\f3\fs18 recombination()
+\f4\fs20 callback. When only
+\f3\fs18 strand1
+\f4\fs20 and
+\f3\fs18 strand3
+\f4\fs20 are supplied, it is very similar to
+\f3\fs18 addCloned()
+\f4\fs20 , creating a clonal offspring, except that the two parental haplosomes need not belong to the same individual (whatever that might mean biologically). Supplying only
+\f3\fs18 strand1
+\f4\fs20 is useful for modeling clonally reproducing haploids, or any chromosome type that is intrinsically haploid, such as the Y chromosome. For a model of clonally reproducing haploids that undergo horizontal gene transfer (HGT), supplying only
+\f3\fs18 strand1
+\f4\fs20 and
+\f3\fs18 strand2
+\f4\fs20 will allow HGT from
+\f3\fs18 strand2
+\f4\fs20 to replace segments of an otherwise clonal copy of
+\f3\fs18 strand1
+\f4\fs20 , while the second haplosome of the generated offspring will be a null haplosome; this could be useful for modeling bacterial conjugation, for example. Other variations are also possible.\
The
\f3\fs18 sex
-\f4\fs20 parameter is interpreted exactly as in
-\f3\fs18 addCrossed()
-\f4\fs20 ; see that method for discussion. If the offspring sex is specified in any way (i.e., if
+\f4\fs20 parameter optionally specifies the sex of the offspring. The default value of
+\f3\fs18 NULL
+\f4\fs20 for
\f3\fs18 sex
-\f4\fs20 is non-
+\f4\fs20 specifies \'93default behavior\'94; in a non-sexual model this is the only legal value, and produces a hermaphroditic offspring. In a sexual model, the \'93default behavior\'94 of
\f3\fs18 NULL
-\f4\fs20 ), the strands provided must be compatible with the sex chosen. If the offspring sex is not specified (i.e., if
+\f4\fs20 is that the offspring\'92s sex is dictated by the haplosome structure it inherits. For example, if the supplied strands indicate that the offspring will have two non-null haplosomes for a chromosome of type
+\f3\fs18 "X"
+\f4\fs20 , the sex of the offspring must therefore be female, whereas if it will have one non-null
+\f3\fs18 "X"
+\f4\fs20 haplosome and one null
+\f3\fs18 "X"
+\f4\fs20 haplosome, it must therefore be male. SLiM will examine the supplied strands to determine such constraints and enforce them. The constraints defined by each chromosome type, as described in
+\f3\fs18 initializeChromosomeType()
+\f4\fs20 , must always be followed, even when using
+\f3\fs18 addRecombinant()
+\f4\fs20 . If
\f3\fs18 sex
\f4\fs20 is
\f3\fs18 NULL
-\f4\fs20 ), the sex will be inferred from the strands provided where possible (when modeling a chromosome with a sex-specific inheritance pattern), or will be chosen randomly otherwise (e.g., when modeling autosomes); it will
-\f1\i not
-\f4\i0 be inferred from the sex of the individuals possessing the parental strands, even when the reproductive mode is essentially clonal from a single parent, since such inference would be ambiguous in the general case.\
-By default, the offspring is considered to have no parents for the purposes of pedigree tracking, since there may be more than two \'93parents\'94 in the general case. If pedigree tracking of parentage is desired,
+\f4\fs20 and the sex of the offspring is unconstrained, the offspring will be chosen as male or female with equal probability. A value of
+\f3\fs18 "M"
+\f4\fs20 or
+\f3\fs18 "F"
+\f4\fs20 for
+\f3\fs18 sex
+\f4\fs20 specifies that the offspring should be male or female, respectively. A
+\f3\fs18 float
+\f4\fs20 value from
+\f3\fs18 0.0
+\f4\fs20 to
+\f3\fs18 1.0
+\f4\fs20 for
+\f3\fs18 sex
+\f4\fs20 provides the probability that the offspring will be male; a value of
+\f3\fs18 0.0
+\f4\fs20 will produce a female, a value of
+\f3\fs18 1.0
+\f4\fs20 will produce a male, and for intermediate values SLiM will draw the sex of the offspring randomly according to the specified probability. In these cases where
+\f3\fs18 sex
+\f4\fs20 is not
+\f3\fs18 NULL
+\f4\fs20 , SLiM will first determine the sex of the individual as just described, and will then examine the supplied strands to confirm that it is compatible with the sex that was determined. Again, if there is a conflict an error will be raised; you cannot specify the sex of an individual to be incompatible with the haplosomes that it inherits, and if you specify a sex with a
+\f3\fs18 string
+\f4\fs20 or
+\f3\fs18 float
+\f4\fs20 value it is up to you to ensure that that is compatible with the supplied strands. (If you need more flexibility, you should probably not use a sexual model at all, but simply use chromosome type
+\f3\fs18 "A"
+\f4\fs20 or
+\f3\fs18 "H"
+\f4\fs20 in a non-sexual model, track the sex of individuals yourself with a tag value such as
+\f3\fs18 tagL0
+\f4\fs20 , and manipulate haplosomes during reproduction however you wish; SLiM then imposes no constraints.)\
+By default, the offspring is considered to have no parents, since there may be more than two \'93parents\'94 in the general case. If specifying parentage is desired,
\f3\fs18 parent1
\f4\fs20 and/or
\f3\fs18 parent2
-\f4\fs20 may be passed to explicitly establish particular individuals as the parents of the offspring for purposes of pedigree tracking. In this case, if only one of
+\f4\fs20 may be passed to explicitly; this will establish those individuals as the parents of the offspring for purposes of pedigree tracking, and for several other purposes described below. If only one of
\f3\fs18 parent1
\f4\fs20 and
\f3\fs18 parent2
@@ -9959,18 +10099,24 @@ By default, the offspring is considered to have no parents for the purposes of p
\f3\fs18 parent1
\f4\fs20 or
\f3\fs18 parent2
-\f4\fs20 to actually be a genetic parent of the offspring at all, although typically they would be.\
+\f4\fs20 to actually be a genetic parent of the offspring at all, although typically they would be. To benefit from the full functionality of
+\f3\fs18 addRecombinant()
+\f4\fs20 as described below, it is best to supply
+\f3\fs18 parent1
+\f4\fs20 and
+\f3\fs18 parent2
+\f4\fs20 when possible.\
The
\f3\fs18 randomizeStrands
\f4\fs20 parameter is used to control the recombination behavior of
\f3\fs18 addRecombinant()
-\f4\fs20 . Two parental strands with crossover breakpoints between them can be specified as described above to generate one offspring strand with recombination \'96 for example, with the
+\f4\fs20 . As described above, two parental strands with crossover breakpoints between them can be specified to generate one offspring strand with recombination \'96 for example, with the
\f3\fs18 strand1
\f4\fs20 ,
\f3\fs18 strand2
\f4\fs20 , and
\f3\fs18 breaks1
-\f4\fs20 parameters, as described above, but the same is true of the
+\f4\fs20 parameters, but the same is true of the
\f3\fs18 strand3
\f4\fs20 ,
\f3\fs18 strand4
@@ -9980,11 +10126,13 @@ The
\f3\fs18 randomizeStrands
\f4\fs20 is
\f3\fs18 F
-\f4\fs20 , the supplied strands are used precisely as given; for example,
+\f4\fs20 , the supplied strands are used as given; for example,
\f3\fs18 strand1
\f4\fs20 will be the initial copy strand when generating the first gamete to form the offspring. This mode should be used if you want explicit control over the initial copy strand; one example would be if your script is explicitly generating all four of the products of a meiosis event. If
\f3\fs18 randomizeStrands
-\f4\fs20 is T, then if
+\f4\fs20 is
+\f3\fs18 T
+\f4\fs20 , then if
\f3\fs18 strand1
\f4\fs20 and
\f3\fs18 strand2
@@ -10018,67 +10166,11 @@ The
\f3\fs18 F
\f4\fs20 prior to SLiM 5, but
\f3\fs18 T
-\f4\fs20 is almost always the correct behavior, as explained above. To try to prevent accidental bugs, the new policy was thus adopted to force the user to explicitly choose
+\f4\fs20 is almost always the correct behavior, as explained above. To try to prevent accidental bugs, this new policy was adopted to force the user to explicitly choose
\f3\fs18 T
\f4\fs20 or
\f3\fs18 F
-\f4\fs20 .)\
-If
-\f3\fs18 randomizeStrands
-\f4\fs20 is
-\f3\fs18 F
-\f4\fs20 (the default),
-\f3\fs18 strand1
-\f4\fs20 will be the initial copy strand when generating the first gamete to form the offspring, and
-\f3\fs18 strand3
-\f4\fs20 will be the initial copy strand when generating the second gamete. If
-\f3\fs18 randomizeStrands
-\f4\fs20 is T, then if
-\f3\fs18 strand1
-\f4\fs20 and
-\f3\fs18 strand2
-\f4\fs20 are both non-
-\f3\fs18 NULL
-\f4\fs20 , 50% of the time they will be swapped, making
-\f3\fs18 strand2
-\f4\fs20 the initial copy strand for the first gamete; and similarly, if
-\f3\fs18 strand3
-\f4\fs20 and
-\f3\fs18 strand4
-\f4\fs20 are both non-
-\f3\fs18 NULL
-\f4\fs20 , 50% of the time they will be swapped, making
-\f3\fs18 strand4
-\f4\fs20 the initial copy strand for the second gamete. This is probably usually the desired behavior, to avoid an inheritance bias due to a lack of randomization in the initial copy strand, so passing
-\f3\fs18 T
-\f4\fs20 for
-\f3\fs18 randomizeStrands
-\f4\fs20 is recommended unless you specifically desire otherwise. It is not the default behavior only for reasons of backward compatibility.\
-These semantics allow several uses for
-\f3\fs18 addRecombinant()
-\f4\fs20 . When all strands are non-
-\f3\fs18 NULL
-\f4\fs20 , it is similar to
-\f3\fs18 addCrossed()
-\f4\fs20 except that the recombination breakpoints are specified explicitly, allowing very precise offspring generation without having to override SLiM\'92s breakpoint generation with a
-\f3\fs18 recombination()
-\f4\fs20 callback. When only
-\f3\fs18 strand1
-\f4\fs20 and
-\f3\fs18 strand3
-\f4\fs20 are supplied, it is very similar to
-\f3\fs18 addCloned()
-\f4\fs20 , creating a clonal offspring, except that the two parental haplosomes need not belong to the same individual (whatever that might mean biologically). Supplying only
-\f3\fs18 strand1
-\f4\fs20 is useful for modeling clonally reproducing haploids; the second haplosome of every offspring will be kept empty and will not receive new mutations. For a model of clonally reproducing haploids that undergo horizontal gene transfer (HGT), supplying only
-\f3\fs18 strand1
-\f4\fs20 and
-\f3\fs18 strand2
-\f4\fs20 will allow HGT from
-\f3\fs18 strand2
-\f4\fs20 to replace segments of an otherwise clonal copy of
-\f3\fs18 strand1
-\f4\fs20 , while the second haplosome of the generated offspring will again be kept empty; this could be useful for modeling bacterial conjugation, for example. Other variations are also possible.\
+\f4\fs20 whenever it matters.)\
The value of the
\f3\fs18 meanParentAge
\f4\fs20 property of the generated offspring is calculated from the mean parent age of each of its two haplosomes (whether they turn out to be null haplosomes or not); that may be an average of two values (if both offspring haplosomes have at least one parent), a single value (if one offspring haplosome has no parent), or no values (if both offspring haplosomes have no parent, in which case
@@ -10091,10 +10183,47 @@ The value of the
\f3\fs18 NULL
\f4\fs20 then that offspring haplosome is parentless and is not used in the final calculation. In other words, if one offspring haplosome has two parents with ages A and B, and the other offspring haplosome has one parent with age C, the
\f3\fs18 meanParentAge
-\f4\fs20 of the offspring will be (A+B+C+C)\'a0/\'a04, not (A+B+C)\'a0/\'a03.\
-Note that gene conversion tracts are not explicitly supported by this method; the
-\f3\fs18 breaks
-\f4\fs20 vectors provide crossover breakpoints, which may be used to implement crossovers or simple gene conversion tracts. There is no way to specify complex gene conversion tracts with heteroduplex mismatch repair.\
+\f4\fs20 of the offspring will be (A+B+C+C)\'a0/\'a04, or equivalently, ((A+B)/2\'a0+\'a0C)\'a0/\'a02, not (A+B+C)\'a0/\'a03.\
+Callbacks can be involved in offspring generation with
+\f3\fs18 addRecombinant()
+\f4\fs20 , but because there are up to four strands with up to four different parents, things are a bit complicated and different from other
+\f3\fs18 add...()
+\f4\fs20 methods; the policy described here seems like the best compromise. The target subpopulation for the
+\f3\fs18 addRecombinant()
+\f4\fs20 call will be used to locate applicable
+\f3\fs18 mutation()
+\f4\fs20 and
+\f3\fs18 modifyChild()
+\f4\fs20 callbacks governing the generation of the offspring individual. On the other hand,
+\f3\fs18 recombination()
+\f4\fs20 callbacks will be found based upon the subpopulations to which
+\f3\fs18 parent1
+\f4\fs20 and
+\f3\fs18 parent2
+\f4\fs20 belong (for reasons discussed further in
+\f3\fs18 drawBreakpoints()
+\f4\fs20 ); if a parent individual is not supplied,
+\f3\fs18 recombination()
+\f4\fs20 callbacks will not be called at all when generating the corresponding offspring haplosome.\
+When breakpoints are explicitly supplied to
+\f3\fs18 addRecombinant()
+\f4\fs20 with
+\f3\fs18 breaks1
+\f4\fs20 or
+\f3\fs18 breaks2
+\f4\fs20 , gene conversion tracts are not well-supported by this method; the
+\f3\fs18 breaks1
+\f4\fs20 and
+\f3\fs18 breaks2
+\f4\fs20 vectors provide simple crossover breakpoints, which may be used to implement crossovers or simple gene conversion tracts, but complex gene conversion tracts with heteroduplex mismatch repair are not supported in this mode of operation since there is no way to supply the relevant information. If, on the other hand,
+\f3\fs18 breaks1
+\f4\fs20 or
+\f3\fs18 breaks2
+\f4\fs20 is
+\f3\fs18 NULL
+\f4\fs20 when generating a haplosome with recombination, then as described above,
+\f3\fs18 addRecombinant()
+\f4\fs20 will generate breakpoints internally for that cross, and in this case, complex gene conversion tracts with heteroduplex mismatch repair are supported, since all of the necessary information is available.\
Beginning in SLiM 4.1, the
\f3\fs18 count
\f4\fs20 parameter dictates how many offspring will be generated (previously, exactly one offspring was generated). Each offspring is generated independently, based upon the given parameters. The returned vector contains all generated offspring, except those that were rejected by a
@@ -10137,18 +10266,18 @@ Also beginning in SLiM 4.1, in spatial models the spatial position of the offspr
\f3\fs18 parent1
\f4\fs20 is
\f3\fs18 NULL
-\f4\fs20 (the default),
+\f4\fs20 ,
\f3\fs18 parent2
\f4\fs20 will be used; if it is also
\f3\fs18 NULL
\f4\fs20 , no spatial position will be inherited.\
Note that this method is only for use in nonWF models. See
\f3\fs18 addCrossed()
-\f4\fs20 for further general notes on the addition of new offspring individuals.\expnd0\expndtw0\kerning0
-\
+\f4\fs20 for further general notes on the addition of new offspring individuals.\
\pard\pardeftab397\li720\fi-446\ri720\sb180\sa60\partightenfactor0
-\f3\fs18 \cf2 \'96\'a0(object)addSelfed(object$\'a0parent, [integer$\'a0count\'a0=\'a01], [logical$\'a0defer\'a0=\'a0F])\
+\f3\fs18 \cf2 \expnd0\expndtw0\kerning0
+\'96\'a0(object)addSelfed(object$\'a0parent, [integer$\'a0count\'a0=\'a01], [logical$\'a0defer\'a0=\'a0F])\
\pard\pardeftab720\li547\ri720\sb60\sa60\partightenfactor0
\f4\fs20 \cf2 \kerning1\expnd0\expndtw0 Generates a new offspring individual from the given parent by selfing, queues it for addition to the target subpopulation, and returns it. The new offspring will not be visible as a member of the target subpopulation until the end of the offspring generation tick cycle stage. The subpopulation of
diff --git a/VERSIONS b/VERSIONS
index a9ff806e..b67b61e0 100644
--- a/VERSIONS
+++ b/VERSIONS
@@ -93,6 +93,9 @@ development head (in the master branch):
policy change: randomizeStrands must now usually be explicitly specified for both addRecombinant() and addMultiRecombinant();
the default is now NULL, and NULL errors so that T or F must be explicitly given unless it does not matter (i.e., there is no recombination)
implement addMultiRecombinant() for complete control over reproduction in multi-chromosome models
+ extend addRecombinant() and addMultiRecombinant() to allow NULL as a breaks vector with two strands, meaning "draw breakpoints for me", to avoid needing drawBreakpoints() in most cases
+ this new code path also allows gene conversion including heteroduplex mismatch repair and gBGC to work with addRecombinant() and addMultiRecombinant()
+ also change to allow a recombination rate to be set for haploid chromosomes; this is potentially useful with addRecombinant() etc.
version 4.3 (Eidos version 3.3):
diff --git a/core/chromosome.cpp b/core/chromosome.cpp
index efa759d5..65c8c5db 100644
--- a/core/chromosome.cpp
+++ b/core/chromosome.cpp
@@ -621,30 +621,6 @@ void Chromosome::ChooseMutationRunLayout(void)
EIDOS_TERMINATION << "ERROR (Chromosome::ChooseMutationRunLayout): (internal error) math error in mutation run calculations." << EidosTerminate();
}
-bool Chromosome::RequiresZeroRecombination(void) const
-{
- switch (type_)
- {
- // these chromosome types allow recombination
- case ChromosomeType::kA_DiploidAutosome:
- case ChromosomeType::kH_HaploidAutosome:
- case ChromosomeType::kX_XSexChromosome:
- case ChromosomeType::kZ_ZSexChromosome:
- return false;
-
- // these types do not, because they are haploid
- case ChromosomeType::kY_YSexChromosome:
- case ChromosomeType::kW_WSexChromosome:
- case ChromosomeType::kHF_HaploidFemaleInherited:
- case ChromosomeType::kFL_HaploidFemaleLine:
- case ChromosomeType::kHM_HaploidMaleInherited:
- case ChromosomeType::kML_HaploidMaleLine:
- case ChromosomeType::kHNull_HaploidAutosomeWithNull:
- case ChromosomeType::kNullY_YSexChromosomeWithNull:
- return true;
- }
-}
-
bool Chromosome::DefaultsToZeroRecombination(void) const
{
switch (type_)
@@ -1407,7 +1383,7 @@ MutationIndex Chromosome::DrawNewMutationExtended(std::pair &p_crossovers) const
+void Chromosome::_DrawCrossoverBreakpoints(IndividualSex p_parent_sex, const int p_num_breakpoints, std::vector &p_crossovers) const
{
// BEWARE! Chromosome::DrawDSBBreakpoints() below must be altered in parallel with this method!
#if DEBUG
@@ -1515,7 +1491,7 @@ void Chromosome::DrawCrossoverBreakpoints(IndividualSex p_parent_sex, const int
// draw a set of uniqued breakpoints according to the "double-stranded break" model and run them through recombination() callbacks, returning the final usable set
// the information returned here also includes a list of heteroduplex regions where mismatches between the two parental strands will need to be resolved
-void Chromosome::DrawDSBBreakpoints(IndividualSex p_parent_sex, const int p_num_breakpoints, std::vector &p_crossovers, std::vector &p_heteroduplex) const
+void Chromosome::_DrawDSBBreakpoints(IndividualSex p_parent_sex, const int p_num_breakpoints, std::vector &p_crossovers, std::vector &p_heteroduplex) const
{
THREAD_SAFETY_IN_ANY_PARALLEL("Chromosome::DrawDSBBreakpoints(): usage of statics, probably many other issues");
@@ -1706,6 +1682,170 @@ void Chromosome::DrawDSBBreakpoints(IndividualSex p_parent_sex, const int p_num_
}
}
+// This high-level function is the funnel for drawing breakpoints. It delegates down to DrawDSBBreakpoints()
+// or DrawCrossoverBreakpoints(), handles recombination() callbacks, and returns a sorted, uniqued vector.
+// You can supply it with a number of breakpoints to draw, or pass -1 to have it draw the number for you.
+// If the caller can handle complex gene conversion tracts, they should pass a vector for those to be placed
+// in. If not, pass nullptr, and this method will raise if complex gene conversion tracts are in use. For
+// addRecombinant() and addMultiRecombinant(), this method allows the parent to be different from the
+// haplosomes that are supplied; the parent individual is used to look up the haplosomes if they are passed
+// as nullptr. The haplosomes are used only if recombination() callbacks are in effect. The parent
+// is also used to look up the sex, for sex-specific recombination rates.
+void Chromosome::DrawBreakpoints(Individual *p_parent, Haplosome *p_haplosome1, Haplosome *p_haplosome2, int p_num_breakpoints, std::vector &p_crossovers, std::vector *p_heteroduplex, const char *p_caller_name)
+{
+ if (!species_.HasGenetics())
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): in " << p_caller_name << ", recombination breakpoints cannot be drawn for a species with no genetics." << EidosTerminate();
+
+ if (!p_heteroduplex && using_DSB_model_ && (simple_conversion_fraction_ != 1.0))
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): in " << p_caller_name << ", complex gene conversion tracts cannot be active since there is no provision for handling heteroduplex regions." << EidosTerminate();
+
+ // look up parent information; note that if parent is nullptr, we do not run recombination() callbacks!
+ // the parent is a required pseudo-parameter for the recombination() callback
+ IndividualSex parent_sex = IndividualSex::kUnspecified;
+ std::vector recombination_callbacks;
+ Subpopulation *parent_subpop = nullptr;
+
+ if (p_parent)
+ {
+ parent_sex = p_parent->sex_;
+ parent_subpop = p_parent->subpopulation_;
+ recombination_callbacks = species_.CallbackBlocksMatching(community_.Tick(), SLiMEidosBlockType::SLiMEidosRecombinationCallback, -1, -1, parent_subpop->subpopulation_id_);
+
+ // SPECIES CONSISTENCY CHECK
+ if (&p_parent->subpopulation_->species_ != &species_)
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): in " << p_caller_name << ", the parent, if supplied, must belong to the same species as the target chromosome." << EidosTerminate();
+ }
+ else // !p_parent
+ {
+ // In a sexual model with sex-specific recombination maps, we need to know the parent we're
+ // generating breakpoints for; in other situations it is optional, but recombination()
+ // breakpoints will not be called if parent is NULL.
+ if (!single_recombination_map_)
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): in " << p_caller_name << ", a parent must be supplied since sex-specific recombination maps are in use (to determine which map to use, from the sex of the parent)." << EidosTerminate();
+ }
+
+ // look up haplosome information, used only for the recombination() callback
+ if ((!p_haplosome1 && p_haplosome2) || (!p_haplosome2 && p_haplosome1))
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): (internal error) in " << p_caller_name << ", haplosomes must either be supplied or not supplied." << EidosTerminate();
+
+ if (p_haplosome1)
+ {
+ // SPECIES CONSISTENCY CHECK
+ if ((&p_haplosome1->OwningIndividual()->subpopulation_->species_ != &species_) ||
+ (&p_haplosome2->OwningIndividual()->subpopulation_->species_ != &species_))
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): in " << p_caller_name << ", parental haplosomes must belong to the same species as the target chromosome." << EidosTerminate();
+
+ if ((p_haplosome1->chromosome_index_ != index_) || (p_haplosome2->chromosome_index_ != index_))
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): in " << p_caller_name << ", parental haplosomes must belong to the target chromosome." << EidosTerminate();
+
+ // Note that if haplosomes were passed in but p_parent is nullptr, we do NOT attempt to infer a parent
+ // subpopulation in order to get recombination() callbacks from it. In the general case that would
+ // not work, since the two haplosomes might belong to different subpopulations; and if we can't do it
+ // in general, we shouldn't try to do it at all. If you want recombination() callbacks, pass p_parent.
+ }
+ else if (p_parent)
+ {
+ // Get the indices of the haplosomes associated with this chromosome. Note that the first/last indices
+ // might be the same, if this is a haploid chromosome. That is OK here. The user is allowed to set a
+ // recombination rate on a haploid chromosome and generate breakpoints for it; what they do with that
+ // information is up to them. (They might use them in an addRecombinant() or addMultiRecombinant() call,
+ // for example.) In that case, of a haploid chromosome, the same single parent haplosome will be passed
+ // twice to recombination() callbacks; that seems better than not defining one of the pseudo-parameters.
+ int first_haplosome_index = species_.FirstHaplosomeIndices()[index_];
+ int last_haplosome_index = species_.LastHaplosomeIndices()[index_];
+
+ // Note that for calling recombination() callbacks below, we treat the parent's first haplosome as the
+ // initial copy strand. If a distinction needs to be made, pass the haplosomes in to this method.
+ p_haplosome1 = p_parent->haplosomes_[first_haplosome_index];
+ p_haplosome2 = p_parent->haplosomes_[last_haplosome_index];
+ }
+
+ // Draw the number of breakpoints, if it was not supplied
+ if (p_num_breakpoints == -1)
+ p_num_breakpoints = DrawBreakpointCount(parent_sex);
+
+ if ((p_num_breakpoints < 0) || (p_num_breakpoints > 1000000L))
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): in " << p_caller_name << ", the number of recombination breakpoints must be in [0, 1000000]." << EidosTerminate();
+
+#if DEBUG
+ if (p_crossovers.size())
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): (internal error) in " << p_caller_name << ", p_crossovers was not supplied empty." << EidosTerminate();
+#endif
+
+ // draw the breakpoints based on the recombination rate map, and sort and unique the result
+ if (p_num_breakpoints)
+ {
+ if (using_DSB_model_)
+ {
+ if (p_heteroduplex)
+ {
+ // p_heteroduplex is not nullptr, so the caller intends to use it for something
+ _DrawDSBBreakpoints(parent_sex, p_num_breakpoints, p_crossovers, *p_heteroduplex);
+ }
+ else
+ {
+ // p_heteroduplex is nullptr, so we need to pass in our own vector; it is not actually used
+ // in this case anyway, since simple_conversion_fraction_ must be 1.0 (as checked above)
+ std::vector heteroduplex;
+
+ _DrawDSBBreakpoints(parent_sex, p_num_breakpoints, p_crossovers, heteroduplex);
+ }
+ }
+ else
+ {
+ _DrawCrossoverBreakpoints(parent_sex, p_num_breakpoints, p_crossovers);
+ }
+
+ // p_crossovers is sorted and uniqued at this point
+
+ if (recombination_callbacks.size())
+ {
+ // a non-zero number of breakpoints, with recombination callbacks
+ bool breaks_changed = species_.population_.ApplyRecombinationCallbacks(p_parent, p_haplosome1, p_haplosome2, p_crossovers, recombination_callbacks);
+
+ // we only sort/unique if the breakpoints have changed, since they were sorted/uniqued before
+ if (breaks_changed && (p_crossovers.size() > 1))
+ {
+ std::sort(p_crossovers.begin(), p_crossovers.end());
+ p_crossovers.erase(unique(p_crossovers.begin(), p_crossovers.end()), p_crossovers.end());
+ }
+ }
+ }
+ else if (recombination_callbacks.size())
+ {
+ // zero breakpoints from the SLiM core, but we have recombination() callbacks
+ species_.population_.ApplyRecombinationCallbacks(p_parent, p_haplosome1, p_haplosome2, p_crossovers, recombination_callbacks);
+
+ if (p_crossovers.size() > 1)
+ {
+ std::sort(p_crossovers.begin(), p_crossovers.end());
+ p_crossovers.erase(unique(p_crossovers.begin(), p_crossovers.end()), p_crossovers.end());
+ }
+ }
+ else
+ {
+ // no breakpoints, no gene conversion, no recombination() callbacks
+ }
+
+ // values in p_crossovers and p_heteroduplex are returned to the caller
+ // p_crossovers is guaranteed to be sorted and uniqued, which we check here
+ // we also check that no position is less than zero or beyond the chromosome end
+#if DEBUG
+ slim_position_t previous_value = -1;
+ slim_position_t last_chrom_position = last_position_;
+
+ for (slim_position_t value : p_crossovers)
+ {
+ if (value <= previous_value)
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): (internal error) breakpoints vector is not sorted/uniqued." << EidosTerminate();
+ if (value > last_chrom_position)
+ EIDOS_TERMINATION << "ERROR (Chromosome::DrawBreakpoints): (internal error) breakpoints vector goes beyond the chromosome end." << EidosTerminate();
+
+ previous_value = value;
+ }
+#endif
+}
+
size_t Chromosome::MemoryUsageForMutationMaps(void)
{
size_t usage = 0;
@@ -2862,60 +3002,17 @@ EidosValue_SP Chromosome::ExecuteMethod_ancestralNucleotides(EidosGlobalStringID
EidosValue_SP Chromosome::ExecuteMethod_drawBreakpoints(EidosGlobalStringID p_method_id, const std::vector &p_arguments, EidosInterpreter &p_interpreter)
{
#pragma unused (p_method_id, p_arguments, p_interpreter)
- if (!species_.HasGenetics())
- EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_drawBreakpoints): drawBreakpoints() may not be called for a species with no genetics." << EidosTerminate();
-
EidosValue *parent_value = p_arguments[0].get();
EidosValue *n_value = p_arguments[1].get();
- if (using_DSB_model_ && (simple_conversion_fraction_ != 1.0))
- EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_drawBreakpoints): drawBreakpoints() does not allow complex gene conversion tracts to be in use, since there is no provision for handling heteroduplex regions." << EidosTerminate();
-
- // In a sexual model with sex-specific recombination maps, we need to know the parent we're
- // generating breakpoints for; in other situations it is optional, but recombination()
- // breakpoints will not be called if parent is NULL.
Individual *parent = nullptr;
if (parent_value->Type() != EidosValueType::kValueNULL)
parent = (Individual *)parent_value->ObjectElementAtIndex_NOCAST(0, nullptr);
- if (!parent && !single_recombination_map_)
- EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_drawBreakpoints): drawBreakpoints() requires a non-NULL parent parameter in sexual models with sex-specific recombination maps." << EidosTerminate();
+ int num_breakpoints = -1; // means "draw them for us"
- IndividualSex parent_sex = IndividualSex::kUnspecified;
- std::vector recombination_callbacks;
- Subpopulation *parent_subpop = nullptr;
-
- // Note that if parent is nullptr, we ignore recombination() callbacks! This is strange, but necessary and documented.
- if (parent)
- {
- parent_sex = parent->sex_;
- parent_subpop = parent->subpopulation_;
- recombination_callbacks = species_.CallbackBlocksMatching(community_.Tick(), SLiMEidosBlockType::SLiMEidosRecombinationCallback, -1, -1, parent_subpop->subpopulation_id_);
-
- // SPECIES CONSISTENCY CHECK
- if (&parent_subpop->species_ != &this->species_)
- EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_drawBreakpoints): drawBreakpoints() requires that parent, if non-NULL, belongs to the same species as the target chromosome." << EidosTerminate();
- }
-
- // Get the indices of the haplosomes associated with this chromosome. Note that the first/last indices might be the same,
- // if this is a haploid chromosome. That is OK here. The user is allowed to set a recombination rate on a haploid
- // chromosome and generate breakpoints for it; what they do with that information is up to them. (They might use them
- // in an addRecombinant() or addMultiRecombinant() call, for example.) In that case, of a haploid chromosome, the same
- // single parent haplosome will be passed twice to recombination() callbacks; that seems better than not defining one of
- // the pseudo-parameters.
- int first_haplosome_index = species_.FirstHaplosomeIndices()[index_];
- int last_haplosome_index = species_.LastHaplosomeIndices()[index_];
-
- // Much of the breakpoint-drawing code here is taken from Population::HaplosomeCrossed().
- // We don't want to split it out into a separate function because (a) we don't want that
- // overhead for HaplosomeCrossed(), which is a hotspot, and (b) we do things slightly
- // differently here (not generating a past-the-end breakpoint, for example).
- int num_breakpoints;
-
- if (n_value->Type() == EidosValueType::kValueNULL)
- num_breakpoints = DrawBreakpointCount(parent_sex);
- else
+ if (n_value->Type() != EidosValueType::kValueNULL)
{
int64_t n = n_value->IntAtIndex_NOCAST(0, nullptr);
@@ -2926,52 +3023,18 @@ EidosValue_SP Chromosome::ExecuteMethod_drawBreakpoints(EidosGlobalStringID p_me
}
std::vector all_breakpoints;
- std::vector heteroduplex; // never actually used since simple_conversion_fraction_ must be 1.0
- // Note that for calling recombination() callbacks below, we always treat the parent's first haplosome as the initial copy strand. This is
- // documented; it is perhaps a weakness of the API here, but if we randomly chose an initial copy strand it would not be used downstream, so.
+ // Note that for calling recombination() callbacks below, we always treat the parent's first haplosome as
+ // the initial copy strand. This is documented; it is perhaps a weakness of the API here, but if we
+ // randomly chose an initial copy strand it would not be used downstream, so.
+ // FIXME an idea: a new parameter, [l$ randomizeStrands = F], could be added that, if true, would
+ // put a breakpoint at 0 half of the time, regardless of recombination rate, so the initial copy
+ // strand is randomized for anyone using the generated breakpoints. This solves the problem, a
+ // little bit clunkily. The main client of this method is users of addRecombinant(), though, and
+ // it now has its own randomizeStrands flag, so maybe this change is unnecessary?
- // draw the breakpoints based on the recombination rate map, and sort and unique the result
- if (num_breakpoints)
- {
- if (using_DSB_model_)
- DrawDSBBreakpoints(parent_sex, num_breakpoints, all_breakpoints, heteroduplex);
- else
- DrawCrossoverBreakpoints(parent_sex, num_breakpoints, all_breakpoints);
-
- if (parent && recombination_callbacks.size())
- {
- // a non-zero number of breakpoints, with recombination callbacks
- Haplosome *parent_haplosome1 = parent->haplosomes_[first_haplosome_index];
- Haplosome *parent_haplosome2 = parent->haplosomes_[last_haplosome_index];
-
- species_.population_.ApplyRecombinationCallbacks(parent_haplosome1, parent_haplosome2, all_breakpoints, recombination_callbacks);
-
- if (all_breakpoints.size() > 1)
- {
- std::sort(all_breakpoints.begin(), all_breakpoints.end());
- all_breakpoints.erase(unique(all_breakpoints.begin(), all_breakpoints.end()), all_breakpoints.end());
- }
- }
- }
- else if (parent && recombination_callbacks.size())
- {
- // zero breakpoints from the SLiM core, but we have recombination() callbacks
- Haplosome *parent_haplosome1 = parent->haplosomes_[first_haplosome_index];
- Haplosome *parent_haplosome2 = parent->haplosomes_[last_haplosome_index];
-
- species_.population_.ApplyRecombinationCallbacks(parent_haplosome1, parent_haplosome2, all_breakpoints, recombination_callbacks);
-
- if (all_breakpoints.size() > 1)
- {
- std::sort(all_breakpoints.begin(), all_breakpoints.end());
- all_breakpoints.erase(unique(all_breakpoints.begin(), all_breakpoints.end()), all_breakpoints.end());
- }
- }
- else
- {
- // no breakpoints, no gene conversion, no recombination() callbacks
- }
+ DrawBreakpoints(parent, /* p_haplosome1 */ nullptr, /* p_haplosome2 */ nullptr, num_breakpoints,
+ all_breakpoints, /* p_heteroduplex */ nullptr, "drawBreakpoints()");
if (all_breakpoints.size() == 0)
return gStaticEidosValue_Integer_ZeroVec;
@@ -3468,9 +3531,6 @@ EidosValue_SP Chromosome::ExecuteMethod_setRecombinationRate(EidosGlobalStringID
if ((recombination_rate < 0.0) || (recombination_rate > 0.5) || std::isnan(recombination_rate))
EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_setRecombinationRate): setRecombinationRate() rate " << recombination_rate << " out of range; rates must be in [0.0, 0.5]." << EidosTerminate();
- if ((recombination_rate != 0.0) && RequiresZeroRecombination())
- EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_setRecombinationRate): setRecombinationRate() requires a recombination rate of 0.0 for all chromosome types except 'A', 'H', 'X', and 'Z'." << EidosTerminate();
-
// then adopt them
rates.clear();
positions.clear();
@@ -3498,9 +3558,6 @@ EidosValue_SP Chromosome::ExecuteMethod_setRecombinationRate(EidosGlobalStringID
if ((recombination_rate < 0.0) || (recombination_rate > 0.5) || std::isnan(recombination_rate))
EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_setRecombinationRate): setRecombinationRate() rate " << recombination_rate << " out of range; rates must be in [0.0, 0.5]." << EidosTerminate();
-
- if ((recombination_rate != 0.0) && RequiresZeroRecombination())
- EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_setRecombinationRate): setRecombinationRate() requires a recombination rate of 0.0 for all chromosome types except 'A', 'H', 'X', and 'Z'." << EidosTerminate();
}
// The stake here is that the last position in the chromosome is not allowed to change after the chromosome is
diff --git a/core/chromosome.h b/core/chromosome.h
index 37e36313..cb028277 100644
--- a/core/chromosome.h
+++ b/core/chromosome.h
@@ -306,7 +306,6 @@ class Chromosome : public EidosDictionaryRetained
inline bool UsingSingleMutationMap(void) const { return single_mutation_map_; }
inline size_t GenomicElementCount(void) const { return genomic_elements_.size(); }
- bool RequiresZeroRecombination(void) const;
bool DefaultsToZeroRecombination(void) const;
// draw the number of mutations that occur, based on the overall mutation rate
@@ -325,9 +324,12 @@ class Chromosome : public EidosDictionaryRetained
// draw the number of breakpoints that occur, based on the overall recombination rate
int DrawBreakpointCount(IndividualSex p_sex) const;
- // choose a set of recombination breakpoints, based on recomb. intervals, overall recomb. rate, and gene conversion parameters
- void DrawCrossoverBreakpoints(IndividualSex p_parent_sex, const int p_num_breakpoints, std::vector &p_crossovers) const;
- void DrawDSBBreakpoints(IndividualSex p_parent_sex, const int p_num_breakpoints, std::vector &p_crossovers, std::vector &p_heteroduplex) const;
+ // choose a set of recombination breakpoints, based on recomb. intervals, overall recomb. rate, and gene
+ // conversion parameters; DrawBreakpoints() is the high-level funnel that most callers should use, whereas
+ // the low-level functions do not handle recombination() callbacks and other niceties
+ void _DrawCrossoverBreakpoints(IndividualSex p_parent_sex, const int p_num_breakpoints, std::vector &p_crossovers) const;
+ void _DrawDSBBreakpoints(IndividualSex p_parent_sex, const int p_num_breakpoints, std::vector &p_crossovers, std::vector &p_heteroduplex) const;
+ void DrawBreakpoints(Individual *p_parent, Haplosome *p_haplosome1, Haplosome *p_haplosome2, int p_num_breakpoints, std::vector &p_crossovers, std::vector *p_heteroduplex, const char *p_caller_name);
#ifndef USE_GSL_POISSON
// draw both the mutation count and breakpoint count, using a single Poisson draw for speed
diff --git a/core/population.cpp b/core/population.cpp
index 31e13fe6..fe112086 100644
--- a/core/population.cpp
+++ b/core/population.cpp
@@ -2459,7 +2459,7 @@ void Population::EvolveSubpopulation(Subpopulation &p_subpop, bool p_mate_choice
}
// apply recombination() callbacks to a generated child; a return of true means breakpoints were changed
-bool Population::ApplyRecombinationCallbacks(Haplosome *p_haplosome1, Haplosome *p_haplosome2, std::vector &p_crossovers, std::vector &p_recombination_callbacks)
+bool Population::ApplyRecombinationCallbacks(Individual *p_parent, Haplosome *p_haplosome1, Haplosome *p_haplosome2, std::vector &p_crossovers, std::vector &p_recombination_callbacks)
{
THREAD_SAFETY_IN_ANY_PARALLEL("Population::ApplyRecombinationCallbacks(): running Eidos callback");
@@ -2518,16 +2518,13 @@ bool Population::ApplyRecombinationCallbacks(Haplosome *p_haplosome1, Haplosome
// the value objects, and we know that the values we are setting here will not change (the objects
// referred to by the values may change, but the values themselves will not change).
if (recombination_callback->contains_individual_)
- {
- Individual *individual = p_haplosome1->OwningIndividual();
- callback_symbols.InitializeConstantSymbolEntry(gID_individual, individual->CachedEidosValue());
- }
+ callback_symbols.InitializeConstantSymbolEntry(gID_individual, p_parent->CachedEidosValue());
if (recombination_callback->contains_haplosome1_)
callback_symbols.InitializeConstantSymbolEntry(gID_haplosome1, p_haplosome1->CachedEidosValue());
if (recombination_callback->contains_haplosome2_)
callback_symbols.InitializeConstantSymbolEntry(gID_haplosome2, p_haplosome2->CachedEidosValue());
if (recombination_callback->contains_subpop_)
- callback_symbols.InitializeConstantSymbolEntry(gID_subpop, p_haplosome1->OwningIndividual()->subpopulation_->SymbolTableEntry().second);
+ callback_symbols.InitializeConstantSymbolEntry(gID_subpop, p_parent->subpopulation_->SymbolTableEntry().second);
// All the variable entries for the crossovers and gene conversion start/end points
if (recombination_callback->contains_breakpoints_)
@@ -2649,8 +2646,9 @@ void Population::HaplosomeCrossed(Chromosome &p_chromosome, Haplosome &p_child_h
// some behaviors -- which callbacks to use, which recombination/mutation rate to use, subpop of origin for mutations, etc. --
// depend upon characteristics of the first parent, so we fetch the necessary properties here
- Subpopulation *source_subpop = parent_haplosome_1->individual_->subpopulation_;
- IndividualSex parent_sex = parent_haplosome_1->individual_->sex_;
+ Individual *parent_individual = parent_haplosome_1->individual_;
+ Subpopulation *source_subpop = parent_individual->subpopulation_;
+ IndividualSex parent_sex = parent_individual->sex_;
// determine how many mutations and breakpoints we have
int num_mutations, num_breakpoints;
@@ -2682,12 +2680,15 @@ void Population::HaplosomeCrossed(Chromosome &p_chromosome, Haplosome &p_child_h
//std::cout << num_mutations << " mutations, " << num_breakpoints << " breakpoints" << std::endl;
// draw the breakpoints based on the recombination rate map, and sort and unique the result
+ // we don't use Chromosome::DrawBreakpoints(), for speed, but this code mirrors it
if (num_breakpoints)
{
if (p_chromosome.using_DSB_model_)
- p_chromosome.DrawDSBBreakpoints(parent_sex, num_breakpoints, all_breakpoints, heteroduplex);
+ p_chromosome._DrawDSBBreakpoints(parent_sex, num_breakpoints, all_breakpoints, heteroduplex);
else
- p_chromosome.DrawCrossoverBreakpoints(parent_sex, num_breakpoints, all_breakpoints);
+ p_chromosome._DrawCrossoverBreakpoints(parent_sex, num_breakpoints, all_breakpoints);
+
+ // all_breakpoints is sorted and uniqued at this point
if (f_callbacks && p_recombination_callbacks)
{
@@ -2695,12 +2696,13 @@ void Population::HaplosomeCrossed(Chromosome &p_chromosome, Haplosome &p_child_h
if (p_chromosome.using_DSB_model_ && (p_chromosome.simple_conversion_fraction_ != 1.0))
EIDOS_TERMINATION << "ERROR (Population::HaplosomeCrossed): recombination() callbacks may not be used when complex gene conversion tracts are in use, since recombination() callbacks have no support for heteroduplex regions." << EidosTerminate();
- ApplyRecombinationCallbacks(parent_haplosome_1, parent_haplosome_2, all_breakpoints, *p_recombination_callbacks);
+ bool breaks_changed = ApplyRecombinationCallbacks(parent_individual, parent_haplosome_1, parent_haplosome_2, all_breakpoints, *p_recombination_callbacks);
num_breakpoints = (int)all_breakpoints.size();
if (num_breakpoints)
{
- if (all_breakpoints.size() > 1)
+ // we only sort/unique if the breakpoints have changed, since they were sorted/uniqued before
+ if (breaks_changed && (all_breakpoints.size() > 1))
{
std::sort(all_breakpoints.begin(), all_breakpoints.end());
all_breakpoints.erase(unique(all_breakpoints.begin(), all_breakpoints.end()), all_breakpoints.end());
@@ -2728,7 +2730,7 @@ void Population::HaplosomeCrossed(Chromosome &p_chromosome, Haplosome &p_child_h
if (p_chromosome.using_DSB_model_ && (p_chromosome.simple_conversion_fraction_ != 1.0))
EIDOS_TERMINATION << "ERROR (Population::HaplosomeCrossed): recombination() callbacks may not be used when complex gene conversion tracts are in use, since recombination() callbacks have no support for heteroduplex regions." << EidosTerminate();
- ApplyRecombinationCallbacks(parent_haplosome_1, parent_haplosome_2, all_breakpoints, *p_recombination_callbacks);
+ ApplyRecombinationCallbacks(parent_individual, parent_haplosome_1, parent_haplosome_2, all_breakpoints, *p_recombination_callbacks);
num_breakpoints = (int)all_breakpoints.size();
if (num_breakpoints)
@@ -2757,6 +2759,15 @@ void Population::HaplosomeCrossed(Chromosome &p_chromosome, Haplosome &p_child_h
}
}
+ // A leading zero in the breakpoints vector switches copy strands before copying begins.
+ // In HaplosomeRecombined() we explicitly remove that, to prevent it from being recorded
+ // in the tree sequence. That's a bit trickier to do here, since we've deferred adding the
+ // past-the-end marker above in some cases. It is not essential to do, I think, just a
+ // nicety that makes for a cleaner recorded tree sequence. Also, a leading zero should
+ // be much less common in HaplosomeCrossed(), since the user is not playing games the way
+ // they do with HaplosomeRecombined(). So I'm going to skip it for now, which I think is
+ // basically harmless; but I'll revisit this later. FIXME MULTICHROM
+
// TREE SEQUENCE RECORDING
bool recording_tree_sequence = f_treeseq;
bool recording_tree_sequence_mutations = f_treeseq && species_.RecordingTreeSequenceMutations();
@@ -3803,8 +3814,24 @@ void Population::HaplosomeRecombined(Chromosome &p_chromosome, Haplosome &p_chil
// determine how many mutations we have
int num_mutations = p_chromosome.DrawMutationCount(parent_sex);
- // we need a defined end breakpoint, so we add it now
- p_breakpoints.emplace_back(p_chromosome.last_position_mutrun_ + 1);
+ // A leading zero in the breakpoints vector switches copy strands before copying begins.
+ // We want to handle that up front, primarily because we don't want to record it in treeseq.
+ // FIXME MULTICHROM erase() here can be very expensive; would be good to switch to using a start pointer
+ // and count, which would allow the first element to be removed in O(1) time by advancing the pointer
+ // I've put off doing that change because it echos across many calls and will make big diffs
+ // As a part of this work, RecordNewHaplosome() should get a derived version for the no-recombination case,
+ // especially the null haplosome case; recording those should be fast to avoid bogging down the model
+ while (p_breakpoints.front() == 0)
+ {
+ p_breakpoints.erase(p_breakpoints.begin());
+ std::swap(parent_haplosome_1, parent_haplosome_2);
+ }
+
+ // we need a defined end breakpoint, so we add it now; we don't want more than one, though
+ // in some cases the caller might already have this breakpoint added, so we need to check
+ // (that happens with multiple children in addRecombinant(), for example)
+ if (p_breakpoints.back() <= p_chromosome.last_position_mutrun_)
+ p_breakpoints.emplace_back(p_chromosome.last_position_mutrun_ + 1);
// TREE SEQUENCE RECORDING
bool recording_tree_sequence = f_treeseq;
diff --git a/core/population.h b/core/population.h
index 4bd9d9f0..ffebf8b9 100644
--- a/core/population.h
+++ b/core/population.h
@@ -220,7 +220,8 @@ class Population
bool ApplyModifyChildCallbacks(Individual *p_child, Individual *p_parent1, Individual *p_parent2, bool p_is_selfing, bool p_is_cloning, Subpopulation *p_target_subpop, Subpopulation *p_source_subpop, std::vector &p_modify_child_callbacks);
// apply recombination() callbacks to a generated child; a return of true means the breakpoints were changed
- bool ApplyRecombinationCallbacks(Haplosome *p_haplosome1, Haplosome *p_haplosome2, std::vector &p_crossovers, std::vector &p_recombination_callbacks);
+ // note that in the general case (e.g., addRecombinant()), the haplosomes might not belong to the parent
+ bool ApplyRecombinationCallbacks(Individual *p_parent, Haplosome *p_haplosome1, Haplosome *p_haplosome2, std::vector &p_crossovers, std::vector &p_recombination_callbacks);
// generate a child haplosome from parental haplosome(s), very directly -- no null haplosomes etc., just cross/clone/recombine
// these methods are templated with variants for speed; see also MungeIndividualCrossed() etc.
diff --git a/core/species.cpp b/core/species.cpp
index fa9d9622..24d2ce2d 100644
--- a/core/species.cpp
+++ b/core/species.cpp
@@ -5111,10 +5111,20 @@ void Species::RecordNewHaplosome(std::vector *p_breakpoints, Ha
right = (*p_breakpoints)[i];
tsk_id_t parent = (tsk_id_t) (polarity ? haplosome1TSKID : haplosome2TSKID);
+ polarity = !polarity;
+
+ // Sometimes the user might add a breakpoint at 0, to flip the initial copy strand, as in the meiotic
+ // drive recipe. If they do that, a left==right breakpoint might make it in to here. That would be
+ // a bug in the caller. This has never been seen in the wild, so I'll make it DEBUG only. In non-
+ // DEBUG runs the tree sequence will fail to pass integrity checks, with TSK_ERR_BAD_EDGE_INTERVAL.
+#if DEBUG
+ if (left >= right)
+ EIDOS_TERMINATION << "ERROR (Species::RecordNewHaplosome): (internal error) a left==right breakpoint was passed to RecordNewHaplosome()." << EidosTerminate();
+#endif
+
int ret = tsk_edge_table_add_row(&tsinfo.tables_.edges, left, right, parent, offspringTSKID, NULL, 0);
if (ret < 0) handle_error("tsk_edge_table_add_row", ret);
- polarity = !polarity;
left = right;
}
diff --git a/core/species_eidos.cpp b/core/species_eidos.cpp
index 5d6f6935..589446e4 100644
--- a/core/species_eidos.cpp
+++ b/core/species_eidos.cpp
@@ -667,13 +667,11 @@ EidosValue_SP Species::ExecuteContextFunction_initializeRecombinationRate(const
double recombination_rate = rates_value->NumericAtIndex_NOCAST(0, nullptr);
- // check values
+ // check values; I thought about requiring a rate of 0.0 for all haploid chromosome types, but maybe
+ // the user wants to recombine them sometimes with addRecombinant(), no need to prevent them
if ((recombination_rate < 0.0) || (recombination_rate > 0.5) || std::isnan(recombination_rate))
EIDOS_TERMINATION << "ERROR (Species::ExecuteContextFunction_initializeRecombinationRate): initializeRecombinationRate() requires rates to be in [0.0, 0.5] (" << EidosStringForFloat(recombination_rate) << " supplied)." << EidosTerminate();
- if ((recombination_rate != 0.0) && chromosome->RequiresZeroRecombination())
- EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_setRecombinationRate): setRecombinationRate() requires a recombination rate of 0.0 for all chromosome types except 'A', 'H', 'X', and 'Z'." << EidosTerminate();
-
// then adopt them
rates.clear();
positions.clear();
@@ -701,9 +699,6 @@ EidosValue_SP Species::ExecuteContextFunction_initializeRecombinationRate(const
if ((recombination_rate < 0.0) || (recombination_rate > 0.5) || std::isnan(recombination_rate))
EIDOS_TERMINATION << "ERROR (Species::ExecuteContextFunction_initializeRecombinationRate): initializeRecombinationRate() requires rates to be in [0.0, 0.5] (" << EidosStringForFloat(recombination_rate) << " supplied)." << EidosTerminate();
- if ((recombination_rate != 0.0) && chromosome->RequiresZeroRecombination())
- EIDOS_TERMINATION << "ERROR (Chromosome::ExecuteMethod_setRecombinationRate): setRecombinationRate() requires a recombination rate of 0.0 for all chromosome types except 'A', 'H', 'X', and 'Z'." << EidosTerminate();
-
if (chromosome->extent_immutable_)
{
if ((recombination_end_position <= 0) || (recombination_end_position > chromosome->last_position_))
diff --git a/core/subpopulation.cpp b/core/subpopulation.cpp
index 7198b623..8538378b 100644
--- a/core/subpopulation.cpp
+++ b/core/subpopulation.cpp
@@ -6940,8 +6940,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
if ((strand1 && strand1->IsNull()) || (strand2 && strand2->IsNull()) || (strand3 && strand3->IsNull()) || (strand4 && strand4->IsNull()))
EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): a parental strand for addMultiRecombinant() is a null haplosome, which is not allowed; pass NULL instead." << EidosTerminate();
- // The parental strands must be visible in the subpopulation, and we need to be able to find them
- // to check their sex
+ // The parental strands must be visible in the subpopulation
Individual *strand1_parent = (strand1 ? strand1->individual_ : nullptr);
Individual *strand2_parent = (strand2 ? strand2->individual_ : nullptr);
Individual *strand3_parent = (strand3 ? strand3->individual_ : nullptr);
@@ -6976,7 +6975,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
// Determine whether we are randomizing strands; note this sets randomizeStrands once for each
// iteration, but always to the same value. We need to check each inheritance dictionary against
- // the specified value of randomizeStrands to validate it.
+ // the specified value of randomizeStrands to validate it, so we need it in pass 1.
if (randomizeStrands_value->Type() == EidosValueType::kValueLogical)
{
randomizeStrands = randomizeStrands_value->LogicalData()[0];
@@ -6997,6 +6996,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
bool haplosome2_null = (!strand3 && !strand4);
// Determine whether we're going to make a second haplosome -- if the chromosome type is diploid.
+ // The second haplosome might be a null haplosome; the question is whether to make one at all.
// _ValidateHaplosomesAndChooseSex() does this check, but only in DEBUG; this is our responsibility.
ChromosomeType inheritance_chromosome_type = inheritance_chromosome->Type();
bool make_second_haplosome = false;
@@ -7008,7 +7008,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
make_second_haplosome = true;
if (!haplosome2_null && !make_second_haplosome)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): for chromosome type '" << inheritance_chromosome_type <<"', addMultiRecombinant() requires that the second offspring haplosome is configured to be a null haplosome (since chromosome type '" << inheritance_chromosome_type << "' is haploid)." << EidosTerminate();
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): for chromosome type '" << inheritance_chromosome_type <<"', addMultiRecombinant() requires that the second offspring haplosome is configured to be a null haplosome (since chromosome type '" << inheritance_chromosome_type << "' is intrinsically haploid)." << EidosTerminate();
// If we're generating any null haplosomes, we need to remember that in the Subpopulation state,
// to turn off optimizations. If the chromosome is haploid, we chack only haplosome1_null.
@@ -7018,6 +7018,9 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
// Check that the breakpoint vectors make sense; breakpoints may not be supplied for a NULL pair or
// a half-NULL pair, but must be supplied for a non-NULL pair. BCH 9/20/2021: Added logic here in
// support of the new semantics that (NULL, NULL, NULL) makes a null haplosome, not an empty haplosome.
+ // BCH 1/8/2025: Changed logic here in support of the new semantics that (haplosome, haplosome, NULL)
+ // requests that SLiM generate breakpoints for a cross in the usual way, so the user doesn't have to.
+ //
// First we need to check the types of the supplied breakpoint vectors; again, Eidos has done no
// validation for us! We do not sort/unique/bounds-check the breakpoint vectors here, since we don't
// need them until the second pass.
@@ -7034,20 +7037,20 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
int breaks1count = breaks1_value->Count(), breaks2count = breaks2_value->Count();
- if (haplosome1_null && (breaks1count != 0))
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): with a NULL strand1 and strand2, breaks1 must be NULL or empty." << EidosTerminate();
- else if ((breaks1count != 0) && !strand2)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): non-empty breaks1 supplied with a NULL strand2; recombination between strand1 and strand2 is not possible, so breaks1 must be NULL or empty." << EidosTerminate();
-
- if (haplosome2_null && (breaks2count != 0))
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): with a NULL strand3 and strand4, breaks2 must be NULL or empty." << EidosTerminate();
- else if ((breaks2count != 0) && !strand4)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): non-empty breaks2 supplied with a NULL strand4; recombination between strand3 and strand4 is not possible, so breaks2 must be NULL or empty." << EidosTerminate();
-
- if ((breaks1_value->Type() == EidosValueType::kValueNULL) && strand1 && strand2)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): strand1 and strand2 are both supplied, so breaks1 may not be NULL (but may be empty)." << EidosTerminate();
- if ((breaks2_value->Type() == EidosValueType::kValueNULL) && strand3 && strand4)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): strand3 and strand4 are both supplied, so breaks2 may not be NULL (but may be empty)." << EidosTerminate();
+ if (breaks1count != 0)
+ {
+ if (haplosome1_null)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): with a NULL strand1 and strand2, breaks1 must be NULL or empty." << EidosTerminate();
+ else if (!strand2)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): non-empty breaks1 supplied with a NULL strand2; recombination between strand1 and strand2 is not possible, so breaks1 must be NULL or empty." << EidosTerminate();
+ }
+ if (breaks2count != 0)
+ {
+ if (haplosome2_null)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): with a NULL strand3 and strand4, breaks2 must be NULL or empty." << EidosTerminate();
+ else if (!strand4)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): non-empty breaks2 supplied with a NULL strand4; recombination between strand3 and strand4 is not possible, so breaks2 must be NULL or empty." << EidosTerminate();
+ }
// The mean parent age is averaged across the mean parent age for each non-null child haplosome
if (strand1_parent && strand2_parent)
@@ -7101,6 +7104,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
if (constrained_child_sex == IndividualSex::kUnspecified)
{
// So far we have seen no constraint; accept the choice _ValidateHaplosomesAndChooseSex() made.
+ // Note that this might still be IndividualSex::kUnspecified, until we get to a constraint.
constrained_child_sex = inheritance_child_sex;
}
else if ((inheritance_child_sex != IndividualSex::kUnspecified) && (constrained_child_sex != inheritance_child_sex))
@@ -7115,8 +7119,10 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
// Figure out the parents for purposes of pedigree recording. If only one parent was supplied, use it
// for both, just as we do for cloning and selfing; it makes relatedness() work. Note mean_parent_age
- // comes from the strands. BCH 9/26/2023 the first parent can now also be used for spatial positioning,
- // even if pedigree tracking is not enabled.
+ // comes from the strands. BCH 9/26/2023: The first parent can now also be used for spatial positioning,
+ // even if pedigree tracking is not enabled. BCH 1/10/2025: We infer selfing/cloning for the parental
+ // pattern here, which addRecombinant() does not do, because we need that to infer missing inheritance
+ // dictionaries below.
bool pedigrees_enabled = species_.PedigreesEnabled();
Individual *pedigree_parent1 = nullptr;
Individual *pedigree_parent2 = nullptr;
@@ -7127,12 +7133,6 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
if (parent2_value->Type() != EidosValueType::kValueNULL)
pedigree_parent2 = (Individual *)parent2_value->ObjectData()[0];
- if ((&pedigree_parent1->subpopulation_->species_ != &species_) || (&pedigree_parent1->subpopulation_->species_ != &species_))
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): addMultiRecombinant() requires that both parents belong to the same species as the target subpopulation." << EidosTerminate();
-
- if ((pedigree_parent1->index_ == -1) || (pedigree_parent1->index_ == -1))
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): parent1 and parent2 must be visible in a subpopulation (i.e., may not be new juveniles)." << EidosTerminate();
-
if (pedigree_parent1 == pedigree_parent2)
is_selfing = true;
@@ -7147,7 +7147,17 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
is_cloning = true;
}
- // Generate the number of children requested
+ if (pedigree_parent1)
+ {
+ // if we have pedigree parents, we need to sanity-check them
+ if ((&pedigree_parent1->subpopulation_->species_ != &species_) || (&pedigree_parent2->subpopulation_->species_ != &species_))
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): addMultiRecombinant() requires that both parents belong to the same species as the target subpopulation." << EidosTerminate();
+
+ if ((pedigree_parent1->index_ == -1) || (pedigree_parent2->index_ == -1))
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): parent1 and parent2 must be visible in a subpopulation (i.e., may not be new juveniles)." << EidosTerminate();
+ }
+
+ // Generate the number of children requested, using mutation() callbacks from the target subpopulation (this)
Eidos_RNG_State *rng_state = EIDOS_STATE_RNG(omp_get_thread_num());
std::vector *mutation_callbacks = ®istered_mutation_callbacks_;
@@ -7159,13 +7169,13 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
for (int64_t child_index = 0; child_index < child_count; ++child_index)
{
- // Determine the sex of the child if it remains unconstrained
+ // If the child's sex is unconstrained, each child generated draws its sex independently
IndividualSex child_sex = constrained_child_sex;
if (child_sex == IndividualSex::kUnspecified)
child_sex = (Eidos_RandomBool(rng_state) ? IndividualSex::kMale : IndividualSex::kFemale);
- // Make the new individual as a candidate
+ // Make the new individual as a candidate; we make its haplosomes in the pass 2 loop below
Individual *individual = NewSubpopIndividual(/* index */ -1, child_sex, /* age */ 0, /* fitness */ NAN, mean_parent_age);
slim_pedigreeid_t pid = 0;
@@ -7200,6 +7210,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
for (Chromosome *chromosome : chromosomes)
{
+ ChromosomeType chromosome_type = chromosome->Type();
EidosDictionaryUnretained *inheritance = nullptr;
if (pattern_integerKeys)
@@ -7218,7 +7229,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
}
if (!inheritance && !pedigree_parent1)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): addMultiRecombinant() has insufficient information to handle a chromosome (id " << chromosome->ID() << ", symbol '" << chromosome->Symbol() << "'); no inheritance dictionary was specified for this chromosome, and parent1 / parent2 are NULL." << EidosTerminate();
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): addMultiRecombinant() has insufficient information to handle a chromosome (id " << chromosome->ID() << ", symbol '" << chromosome->Symbol() << "'); no inheritance dictionary was specified for this chromosome, and inference of a default behavior is not possible since parent1 and parent2 are NULL." << EidosTerminate();
// this is the information we will need, from the inheritance dictionary or synthesized
Haplosome *strand1 = nullptr;
@@ -7226,42 +7237,50 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
Haplosome *strand3 = nullptr;
Haplosome *strand4 = nullptr;
std::vector breakvec1, breakvec2;
+ bool generate_breakvec1 = false;
+ bool generate_breakvec2 = false;
if (inheritance)
{
// process the inheritance dictionary; note that validation was done in pass 1
- const EidosDictionaryHashTable_StringKeys *inheritance_stringKeys = inheritance->DictionarySymbols_StringKeys();
- EidosValue *strand1_value = nullptr;
- EidosValue *strand2_value = nullptr;
EidosValue *breaks1_value = nullptr;
- EidosValue *strand3_value = nullptr;
- EidosValue *strand4_value = nullptr;
EidosValue *breaks2_value = nullptr;
- for (auto const &inheritance_element : *inheritance_stringKeys)
{
- const std::string &inheritance_key = inheritance_element.first;
- EidosValue * const inheritance_value = inheritance_element.second.get();
+ const EidosDictionaryHashTable_StringKeys *inheritance_stringKeys = inheritance->DictionarySymbols_StringKeys();
+ EidosValue *strand1_value = nullptr;
+ EidosValue *strand2_value = nullptr;
+ EidosValue *strand3_value = nullptr;
+ EidosValue *strand4_value = nullptr;
- if (inheritance_key == strand1_string) strand1_value = inheritance_value;
- else if (inheritance_key == strand2_string) strand2_value = inheritance_value;
- else if (inheritance_key == breaks1_string) breaks1_value = inheritance_value;
- else if (inheritance_key == strand3_string) strand3_value = inheritance_value;
- else if (inheritance_key == strand4_string) strand4_value = inheritance_value;
- else if (inheritance_key == breaks2_string) breaks2_value = inheritance_value;
- }
-
- // Get the haplosomes for the supplied strands, or nullptr for NULL
- if (strand1_value && (strand1_value->Type() == EidosValueType::kValueObject))
- strand1 = (Haplosome *)strand1_value->ObjectData()[0];
- if (strand2_value && (strand2_value->Type() == EidosValueType::kValueObject))
- strand2 = (Haplosome *)strand2_value->ObjectData()[0];
- if (strand3_value && (strand3_value->Type() == EidosValueType::kValueObject))
- strand3 = (Haplosome *)strand3_value->ObjectData()[0];
- if (strand4_value && (strand4_value->Type() == EidosValueType::kValueObject))
- strand4 = (Haplosome *)strand4_value->ObjectData()[0];
-
- // Sort/unique/bounds-check the breakpoint vectors
+ for (auto const &inheritance_element : *inheritance_stringKeys)
+ {
+ const std::string &inheritance_key = inheritance_element.first;
+ EidosValue * const inheritance_value = inheritance_element.second.get();
+
+ if (inheritance_key == strand1_string) strand1_value = inheritance_value;
+ else if (inheritance_key == strand2_string) strand2_value = inheritance_value;
+ else if (inheritance_key == breaks1_string) breaks1_value = inheritance_value;
+ else if (inheritance_key == strand3_string) strand3_value = inheritance_value;
+ else if (inheritance_key == strand4_string) strand4_value = inheritance_value;
+ else if (inheritance_key == breaks2_string) breaks2_value = inheritance_value;
+ }
+
+ // Get the haplosomes for the supplied strands, or nullptr for NULL
+ if (strand1_value && (strand1_value->Type() == EidosValueType::kValueObject))
+ strand1 = (Haplosome *)strand1_value->ObjectData()[0];
+ if (strand2_value && (strand2_value->Type() == EidosValueType::kValueObject))
+ strand2 = (Haplosome *)strand2_value->ObjectData()[0];
+ if (strand3_value && (strand3_value->Type() == EidosValueType::kValueObject))
+ strand3 = (Haplosome *)strand3_value->ObjectData()[0];
+ if (strand4_value && (strand4_value->Type() == EidosValueType::kValueObject))
+ strand4 = (Haplosome *)strand4_value->ObjectData()[0];
+ }
+
+ // Copy the breakpoints into local buffers, since we will sort and modify them below.
+ // If NULL is supplied for breaks, that is a request for the recombination breakpoints
+ // to be generated automatically for a cross. Note that the breaks count could also be
+ // zero because integer(0) was passed; that just means "no breakpoints for this cross".
if (!breaks1_value)
breaks1_value = gStaticEidosValueNULL.get();
if (!breaks2_value)
@@ -7275,21 +7294,10 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
for (int break_index = 0; break_index < breaks1count; ++break_index)
breakvec1.emplace_back(SLiMCastToPositionTypeOrRaise(breaks1_data[break_index]));
-
- std::sort(breakvec1.begin(), breakvec1.end());
- breakvec1.erase(unique(breakvec1.begin(), breakvec1.end()), breakvec1.end());
-
- if (breakvec1.back() > chromosome->last_position_)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): breaks1 contained a value (" << breakvec1.back() << ") that lies beyond the end of the chromosome." << EidosTerminate();
-
- // handle a breakpoint at position 0, which swaps the initial strand;
- // HaplosomeRecombined() does not like this so we do it here
- if (breakvec1.front() == 0)
- {
- breakvec1.erase(breakvec1.begin());
- std::swap(strand1, strand2);
- //std::swap(strand1_value, strand2_value); // not used henceforth
- }
+ }
+ else if ((breaks1_value->Type() == EidosValueType::kValueNULL) && strand1 && strand2)
+ {
+ generate_breakvec1 = true;
}
if (breaks2count)
@@ -7298,21 +7306,10 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
for (int break_index = 0; break_index < breaks2count; ++break_index)
breakvec2.emplace_back(SLiMCastToPositionTypeOrRaise(breaks2_data[break_index]));
-
- std::sort(breakvec2.begin(), breakvec2.end());
- breakvec2.erase(unique(breakvec2.begin(), breakvec2.end()), breakvec2.end());
-
- if (breakvec2.back() > chromosome->last_position_)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): breaks2 contained a value (" << breakvec2.back() << ") that lies beyond the end of the chromosome." << EidosTerminate();
-
- // handle a breakpoint at position 0, which swaps the initial strand;
- // HaplosomeRecombined() does not like this so we do it here
- if (breakvec2.front() == 0)
- {
- breakvec2.erase(breakvec2.begin());
- std::swap(strand3, strand4);
- //std::swap(strand3_value, strand4_value); // not used henceforth
- }
+ }
+ else if ((breaks2_value->Type() == EidosValueType::kValueNULL) && strand3 && strand4)
+ {
+ generate_breakvec2 = true;
}
// Randomly swap initial copy strands, if requested and applicable; this should not alter
@@ -7320,28 +7317,60 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
if (randomizeStrands)
{
if (strand1 && strand2 && Eidos_RandomBool(rng_state))
- {
std::swap(strand1, strand2);
- //std::swap(strand1_value, strand2_value); // not used henceforth
- }
if (strand3 && strand4 && Eidos_RandomBool(rng_state))
- {
std::swap(strand3, strand4);
- //std::swap(strand3_value, strand4_value); // not used henceforth
- }
}
}
else if (is_cloning)
{
- // Infer the inheritance dictionary information as addPatternForClone() would
+ // Infer the inheritance dictionary info as addPatternForClone() would (same code path).
#warning implement me!
+ //species_.InferInheritanceForClone(chromosome_type, pedigree_parent1, child_sex, &strand1, &strand2, breakvec1, &generate_breakvec1);
}
- else // crossed and selfed cases; not sure if selfing needs any special handling here
+ else // crossed and selfed cases
{
- // Infer the inheritance dictionary information as addPatternForCross() would
+ // Infer the inheritance dictionary info as addPatternForCross() would (same code path).
+ // Note that this randomizes the strand order for us, so we can ignore randomizeStrands.
#warning implement me!
+ //species_.InferInheritanceForCross(chromosome_type, pedigree_parent1, pedigree_parent2, child_sex, &strand1, &strand2, breakvec1, &generate_breakvec1, &strand3, &strand4, breakvec2, &generate_breakvec2);
+ }
+
+ // Generate breakpoints. If the user passed NULL for breaks1, but strand1 and strand2 are
+ // both non-NULL, it used to be an error but in SLiM 5 it requests that SLiM generate the
+ // breakpoints for the cross automatically, avoiding the need to call drawBreakpoints().
+ // This also allows addRecombinant() and addMultiRecombinant() to do gene conversion with
+ // heteroduplex mismatch repair and gBGC, which was not previously possible.
+ std::vector heteroduplex1, heteroduplex2;
+
+ if (generate_breakvec1)
+ {
+ chromosome->DrawBreakpoints(pedigree_parent1, strand1, strand2, /* p_num_breakpoints */ -1, breakvec1, &heteroduplex1, "addMultiRecombinant()");
+ }
+ else
+ {
+ std::sort(breakvec1.begin(), breakvec1.end());
+ breakvec1.erase(unique(breakvec1.begin(), breakvec1.end()), breakvec1.end());
+ }
+
+ if (generate_breakvec2)
+ {
+ chromosome->DrawBreakpoints(pedigree_parent2, strand3, strand4, /* p_num_breakpoints */ -1, breakvec2, &heteroduplex2, "addMultiRecombinant()");
+ }
+ else
+ {
+ std::sort(breakvec2.begin(), breakvec2.end());
+ breakvec2.erase(unique(breakvec2.begin(), breakvec2.end()), breakvec2.end());
}
+ if (breakvec1.back() > chromosome->last_position_)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): breaks1 contained a value (" << breakvec1.back() << ") that lies beyond the end of the chromosome." << EidosTerminate();
+ if (breakvec2.back() > chromosome->last_position_)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addMultiRecombinant): breaks2 contained a value (" << breakvec2.back() << ") that lies beyond the end of the chromosome." << EidosTerminate();
+
+ // We used to need to look for a leading 0 in the breaks vectors, and swap the corresponding strands,
+ // but HaplosomeRecombined() now handles that for us since it is shared functionality.
+
//
// Now, one way or another, we have all the information we need to generate the offspring
// haplosomes for the current chromosome. Note that for the selfed/crossed case this still
@@ -7363,7 +7392,6 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
// Determine whether we're going to make a second haplosome -- if the chromosome type is diploid.
// We do this also in pass 1, but we need to do it again because we need make_second_haplosome.
- ChromosomeType chromosome_type = chromosome->Type();
bool make_second_haplosome = false;
if ((chromosome_type == ChromosomeType::kA_DiploidAutosome) ||
@@ -7377,7 +7405,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
Haplosome *haplosome2 = nullptr;
if (make_second_haplosome)
- haplosome2_null ? chromosome->NewHaplosome_NULL(individual, 1) : chromosome->NewHaplosome_NONNULL(individual, 1);
+ haplosome2 = haplosome2_null ? chromosome->NewHaplosome_NULL(individual, 1) : chromosome->NewHaplosome_NONNULL(individual, 1);
if (pedigrees_enabled)
{
@@ -7409,6 +7437,9 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
else
{
(population_.*(population_.HaplosomeRecombined_TEMPLATED))(*chromosome, *haplosome1, strand1, strand2, breakvec1, mutation_callbacks);
+
+ if (heteroduplex1.size() > 0)
+ population_.DoHeteroduplexRepair(heteroduplex1, breakvec1, strand1, strand2, haplosome1);
}
}
else
@@ -7453,6 +7484,9 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
else
{
(population_.*(population_.HaplosomeRecombined_TEMPLATED))(*chromosome, *haplosome2, strand3, strand4, breakvec2, mutation_callbacks);
+
+ if (heteroduplex2.size() > 0)
+ population_.DoHeteroduplexRepair(heteroduplex2, breakvec2, strand3, strand4, haplosome2);
}
}
else
@@ -7469,9 +7503,9 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
}
}
}
- else
+ else if (make_second_haplosome)
{
- // both strands are NULL, so we make a null haplosome; we do nothing but record it
+ // both strands are NULL and this is a diploid chromosome, so we record a null haplosome
if (species_.RecordingTreeSequence())
species_.RecordNewHaplosome(nullptr, haplosome2, nullptr, nullptr);
@@ -7485,8 +7519,6 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
}
// Run the candidate past modifyChild() callbacks; the target subpop's registered callbacks are used
- bool proposed_child_accepted = true;
-
if (registered_modify_child_callbacks_.size())
{
// BCH 4/5/2022: When removing excess pseudo-parameters from callbacks, we lost a bit of
@@ -7495,12 +7527,10 @@ EidosValue_SP Subpopulation::ExecuteMethod_addMultiRecombinant(EidosGlobalString
// documented, and I doubt anybody was using it, and they can do the same without the
// modifyChild() callback, so I'm not viewing this loss of functionality as an obstacle
// to making this change.
- // In addMultiRecombinant() we follow addRecombinant()'s lead: we pass nullptr for the
- // parents, and false for p_is_selfing and p_is_cloning. We could use pedigree_parent1
- // and pedigree_parent2, and draw inferences from them about selfing/cloning, but that
- // would not always be correct; it would just be a guess. The user knows what they are
- // doing, and can script what they need to do.
- proposed_child_accepted = population_.ApplyModifyChildCallbacks(individual, /* p_parent1 */ nullptr, /* p_parent2 */ nullptr, /* p_is_selfing */ false, /* p_is_cloning */ false, /* p_target_subpop */ this, /* p_source_subpop */ nullptr, registered_modify_child_callbacks_);
+ // BCH 1/10/2025: Note that we pass nullptr for the parents, and false for p_is_selfing and
+ // p_is_cloning. We could use pedigree_parent1 and pedigree_parent2, and draw inferences from
+ // them about selfing/cloning, but that would not always be correct; it would just be a guess.
+ bool proposed_child_accepted = population_.ApplyModifyChildCallbacks(individual, /* p_parent1 */ nullptr, /* p_parent2 */ nullptr, /* p_is_selfing */ false, /* p_is_cloning */ false, /* p_target_subpop */ this, /* p_source_subpop */ nullptr, registered_modify_child_callbacks_);
if (!proposed_child_accepted)
{
@@ -7591,11 +7621,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
Chromosome *chromosome = species_.Chromosomes()[0];
// Get arguments and do trivial processing
- EidosValue *strand1_value = p_arguments[0].get();
- EidosValue *strand2_value = p_arguments[1].get();
EidosValue *breaks1_value = p_arguments[2].get();
- EidosValue *strand3_value = p_arguments[3].get();
- EidosValue *strand4_value = p_arguments[4].get();
EidosValue *breaks2_value = p_arguments[5].get();
EidosValue *sex_value = p_arguments[6].get();
EidosValue *parent1_value = p_arguments[7].get();
@@ -7628,12 +7654,26 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
// both NULL and breaks1 is NULL/empty, the offspring haplosome will be a *null* haplosome (not just empty),
// and as before will not receive mutations. That is the way it always should have worked. Again, mutatis
// mutandis, for strand3, strand4, and breaks2. See https://github.com/MesserLab/SLiM/issues/205.
+ // BCH 1/8/2025: A further change for SLiM 5: if strand1 and strand2 are both non-NULL (doing a cross),
+ // breaks1 may now be NULL indicating that addRecombinant() should draw breakpoints in the usual way.
// Get the haplosomes for the supplied strands, or nullptr for NULL
- Haplosome *strand1 = ((strand1_value->Type() == EidosValueType::kValueNULL) ? nullptr : (Haplosome *)strand1_value->ObjectData()[0]);
- Haplosome *strand2 = ((strand2_value->Type() == EidosValueType::kValueNULL) ? nullptr : (Haplosome *)strand2_value->ObjectData()[0]);
- Haplosome *strand3 = ((strand3_value->Type() == EidosValueType::kValueNULL) ? nullptr : (Haplosome *)strand3_value->ObjectData()[0]);
- Haplosome *strand4 = ((strand4_value->Type() == EidosValueType::kValueNULL) ? nullptr : (Haplosome *)strand4_value->ObjectData()[0]);
+ Haplosome *strand1;
+ Haplosome *strand2;
+ Haplosome *strand3;
+ Haplosome *strand4;
+
+ {
+ EidosValue *strand1_value = p_arguments[0].get();
+ EidosValue *strand2_value = p_arguments[1].get();
+ EidosValue *strand3_value = p_arguments[3].get();
+ EidosValue *strand4_value = p_arguments[4].get();
+
+ strand1 = ((strand1_value->Type() == EidosValueType::kValueNULL) ? nullptr : (Haplosome *)strand1_value->ObjectData()[0]);
+ strand2 = ((strand2_value->Type() == EidosValueType::kValueNULL) ? nullptr : (Haplosome *)strand2_value->ObjectData()[0]);
+ strand3 = ((strand3_value->Type() == EidosValueType::kValueNULL) ? nullptr : (Haplosome *)strand3_value->ObjectData()[0]);
+ strand4 = ((strand4_value->Type() == EidosValueType::kValueNULL) ? nullptr : (Haplosome *)strand4_value->ObjectData()[0]);
+ }
// New in SLiM 5, we raise if a null haplosome was passed in; remarkably, this was not checked for
// previously, and could lead to a crash if the user tried to do it! It never makes sense to do it,
@@ -7643,8 +7683,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
if ((strand1 && strand1->IsNull()) || (strand2 && strand2->IsNull()) || (strand3 && strand3->IsNull()) || (strand4 && strand4->IsNull()))
EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): a parental strand for addRecombinant() is a null haplosome, which is not allowed; pass NULL instead." << EidosTerminate();
- // The parental strands must be visible in the subpopulation, and we need to be able to find them
- // to check their sex
+ // The parental strands must be visible in the subpopulation
Individual *strand1_parent = (strand1 ? strand1->individual_ : nullptr);
Individual *strand2_parent = (strand2 ? strand2->individual_ : nullptr);
Individual *strand3_parent = (strand3 ? strand3->individual_ : nullptr);
@@ -7671,7 +7710,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
if (!strand3 && strand4)
EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): if strand3 is NULL, strand4 must also be NULL." << EidosTerminate();
- // Determine whether we are randomizing strands
+ // Determine whether we are randomizing strands; if we are, we will do it for each child generated
bool randomizeStrands = false;
if (randomizeStrands_value->Type() == EidosValueType::kValueLogical)
@@ -7694,6 +7733,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
bool haplosome2_null = (!strand3 && !strand4);
// Determine whether we're going to make a second haplosome -- if the chromosome type is diploid.
+ // The second haplosome might be a null haplosome; the question is whether to make one at all.
// _ValidateHaplosomesAndChooseSex() does this check, but only in DEBUG; this is our responsibility.
ChromosomeType chromosome_type = chromosome->Type();
bool make_second_haplosome = false;
@@ -7705,35 +7745,42 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
make_second_haplosome = true;
if (!haplosome2_null && !make_second_haplosome)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): for chromosome type '" << chromosome_type <<"', addRecombinant() requires that the second offspring haplosome is configured to be a null haplosome (since chromosome type '" << chromosome_type << "' is haploid)." << EidosTerminate();
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): for chromosome type '" << chromosome_type <<"', addRecombinant() requires that the second offspring haplosome is configured to be a null haplosome (since chromosome type '" << chromosome_type << "' is intrinsically haploid)." << EidosTerminate();
// If we're generating any null haplosomes, we need to remember that in the Subpopulation state,
- // to turn off optimizations. If the chromosome is haploid, we chack only haplosome1_null.
+ // to turn off optimizations. If the chromosome is haploid, we check only haplosome1_null.
if (haplosome1_null || (haplosome2_null && make_second_haplosome))
has_null_haplosomes_ = true;
// Check that the breakpoint vectors make sense; breakpoints may not be supplied for a NULL pair or
// a half-NULL pair, but must be supplied for a non-NULL pair. BCH 9/20/2021: Added logic here in
// support of the new semantics that (NULL, NULL, NULL) makes a null haplosome, not an empty haplosome.
+ // BCH 1/8/2025: Changed logic here in support of the new semantics that (haplosome, haplosome, NULL)
+ // requests that SLiM generate breakpoints for a cross in the usual way, so the user doesn't have to.
int breaks1count = breaks1_value->Count(), breaks2count = breaks2_value->Count();
- if (haplosome1_null && (breaks1count != 0))
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): with a NULL strand1 and strand2, breaks1 must be NULL or empty." << EidosTerminate();
- else if ((breaks1count != 0) && !strand2)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): non-empty breaks1 supplied with a NULL strand2; recombination between strand1 and strand2 is not possible, so breaks1 must be NULL or empty." << EidosTerminate();
-
- if (haplosome2_null && (breaks2count != 0))
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): with a NULL strand3 and strand4, breaks2 must be NULL or empty." << EidosTerminate();
- else if ((breaks2count != 0) && !strand4)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): non-empty breaks2 supplied with a NULL strand4; recombination between strand3 and strand4 is not possible, so breaks2 must be NULL or empty." << EidosTerminate();
-
- if ((breaks1_value->Type() == EidosValueType::kValueNULL) && strand1 && strand2)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): strand1 and strand2 are both supplied, so breaks1 may not be NULL (but may be empty)." << EidosTerminate();
- if ((breaks2_value->Type() == EidosValueType::kValueNULL) && strand3 && strand4)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): strand3 and strand4 are both supplied, so breaks2 may not be NULL (but may be empty)." << EidosTerminate();
+ if (breaks1count != 0)
+ {
+ if (haplosome1_null)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): with a NULL strand1 and strand2, breaks1 must be NULL or empty." << EidosTerminate();
+ else if (!strand2)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): non-empty breaks1 supplied with a NULL strand2; recombination between strand1 and strand2 is not possible, so breaks1 must be NULL or empty." << EidosTerminate();
+ }
+ if (breaks2count != 0)
+ {
+ if (haplosome2_null)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): with a NULL strand3 and strand4, breaks2 must be NULL or empty." << EidosTerminate();
+ else if (!strand4)
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): non-empty breaks2 supplied with a NULL strand4; recombination between strand3 and strand4 is not possible, so breaks2 must be NULL or empty." << EidosTerminate();
+ }
- // Sort and unique and bounds-check the breakpoints
+ // Copy the breakpoints into local buffers, since we will sort and modify them below.
+ // If NULL is supplied for breaks, that is a request for the recombination breakpoints
+ // to be generated automatically for a cross. Note that the breaks count could also be
+ // zero because integer(0) was passed; that just means "no breakpoints for this cross".
std::vector breakvec1, breakvec2;
+ bool generate_breakvec1 = false;
+ bool generate_breakvec2 = false;
if (breaks1count)
{
@@ -7741,22 +7788,10 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
for (int break_index = 0; break_index < breaks1count; ++break_index)
breakvec1.emplace_back(SLiMCastToPositionTypeOrRaise(breaks1_data[break_index]));
-
- std::sort(breakvec1.begin(), breakvec1.end());
- breakvec1.erase(unique(breakvec1.begin(), breakvec1.end()), breakvec1.end());
-
- if (breakvec1.back() > chromosome->last_position_)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): breaks1 contained a value (" << breakvec1.back() << ") that lies beyond the end of the chromosome." << EidosTerminate();
-
- // handle a breakpoint at position 0, which swaps the initial strand;
- // HaplosomeRecombined() does not like this so we do it here
- if (breakvec1.front() == 0)
- {
- breakvec1.erase(breakvec1.begin());
- std::swap(strand1, strand2);
- std::swap(strand1_parent, strand2_parent);
- //std::swap(strand1_value, strand2_value); // not used henceforth
- }
+ }
+ else if ((breaks1_value->Type() == EidosValueType::kValueNULL) && strand1 && strand2)
+ {
+ generate_breakvec1 = true;
}
if (breaks2count)
@@ -7765,22 +7800,10 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
for (int break_index = 0; break_index < breaks2count; ++break_index)
breakvec2.emplace_back(SLiMCastToPositionTypeOrRaise(breaks2_data[break_index]));
-
- std::sort(breakvec2.begin(), breakvec2.end());
- breakvec2.erase(unique(breakvec2.begin(), breakvec2.end()), breakvec2.end());
-
- if (breakvec2.back() > chromosome->last_position_)
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): breaks2 contained a value (" << breakvec2.back() << ") that lies beyond the end of the chromosome." << EidosTerminate();
-
- // handle a breakpoint at position 0, which swaps the initial strand;
- // HaplosomeRecombined() does not like this so we do it here
- if (breakvec2.front() == 0)
- {
- breakvec2.erase(breakvec2.begin());
- std::swap(strand3, strand4);
- std::swap(strand3_parent, strand4_parent);
- //std::swap(strand3_value, strand4_value); // not used henceforth
- }
+ }
+ else if ((breaks2_value->Type() == EidosValueType::kValueNULL) && strand3 && strand4)
+ {
+ generate_breakvec2 = true;
}
// The mean parent age is averaged across the mean parent age for each non-null child haplosome
@@ -7832,8 +7855,9 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
// Figure out the parents for purposes of pedigree recording. If only one parent was supplied, use it
// for both, just as we do for cloning and selfing; it makes relatedness() work. Note mean_parent_age
- // comes from the strands. BCH 9/26/2023 the first parent can now also be used for spatial positioning,
- // even if pedigree tracking is not enabled.
+ // comes from the strands. BCH 9/26/2023: the first parent can now also be used for spatial positioning,
+ // even if pedigree tracking is not enabled. BCH 1/10/2025: We do not need to infer selfing/cloning
+ // here, the way addMultiRecombinant() does, since we don't need to infer inheritance patterns.
bool pedigrees_enabled = species_.PedigreesEnabled();
Individual *pedigree_parent1 = nullptr;
Individual *pedigree_parent2 = nullptr;
@@ -7843,17 +7867,21 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
if (parent2_value->Type() != EidosValueType::kValueNULL)
pedigree_parent2 = (Individual *)parent2_value->ObjectData()[0];
- if ((&pedigree_parent1->subpopulation_->species_ != &species_) || (&pedigree_parent1->subpopulation_->species_ != &species_))
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): addRecombinant() requires that both parents belong to the same species as the target subpopulation." << EidosTerminate();
-
- if ((pedigree_parent1->index_ == -1) || (pedigree_parent1->index_ == -1))
- EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): parent1 and parent2 must be visible in a subpopulation (i.e., may not be new juveniles)." << EidosTerminate();
-
if (pedigree_parent1 && !pedigree_parent2)
pedigree_parent2 = pedigree_parent1;
if (pedigree_parent2 && !pedigree_parent1)
pedigree_parent1 = pedigree_parent2;
+ if (pedigree_parent1)
+ {
+ // if we have pedigree parents, we need to sanity-check them
+ if ((&pedigree_parent1->subpopulation_->species_ != &species_) || (&pedigree_parent2->subpopulation_->species_ != &species_))
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): addRecombinant() requires that both parents belong to the same species as the target subpopulation." << EidosTerminate();
+
+ if ((pedigree_parent1->index_ == -1) || (pedigree_parent2->index_ == -1))
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): parent1 and parent2 must be visible in a subpopulation (i.e., may not be new juveniles)." << EidosTerminate();
+ }
+
// Figure out what sex the offspring has to be, based on the chromosome type and haplosomes provided.
// If sex_value specifies a sex, or a probability of a sex, then the sex will be chosen and checked
// against the user-supplied data. If it is NULL, then if the data constrains the choice it will be
@@ -7862,7 +7890,7 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
// chromosome type must be followed, because it is assumed in many places.
IndividualSex constrained_child_sex = _ValidateHaplosomesAndChooseSex(chromosome_type, haplosome1_null, haplosome2_null, sex_value, "addRecombinant()");
- // Generate the number of children requested
+ // Generate the number of children requested, using mutation() callbacks from the target subpopulation (this)
Eidos_RNG_State *rng_state = EIDOS_STATE_RNG(omp_get_thread_num());
std::vector *mutation_callbacks = ®istered_mutation_callbacks_;
@@ -7881,30 +7909,68 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
child_sex = (Eidos_RandomBool(rng_state) ? IndividualSex::kMale : IndividualSex::kFemale);
// Randomly swap initial copy strands, if requested and applicable; this should not alter
- // any of the decisions we made earlier about null vs. non-null, child sex, etc.
+ // any of the decisions we made earlier about null vs. non-null, child sex, etc. Note that
+ // with child_count > 1 we'll be swapping back and forth, but that's fine.
if (randomizeStrands)
{
if (strand1 && strand2 && Eidos_RandomBool(rng_state))
{
std::swap(strand1, strand2);
std::swap(strand1_parent, strand2_parent);
- //std::swap(strand1_value, strand2_value); // not used henceforth
}
if (strand3 && strand4 && Eidos_RandomBool(rng_state))
{
std::swap(strand3, strand4);
std::swap(strand3_parent, strand4_parent);
- //std::swap(strand3_value, strand4_value); // not used henceforth
}
}
+ // Generate breakpoints. If the user passed NULL for breaks1, but strand1 and strand2 are
+ // both non-NULL, it used to be an error but in SLiM 5 it requests that SLiM generate the
+ // breakpoints for the cross automatically, avoiding the need to call drawBreakpoints().
+ // This also allows addRecombinant() and addMultiRecombinant() to do gene conversion with
+ // heteroduplex mismatch repair and gBGC, which was not previously possible.
+ std::vector heteroduplex1, heteroduplex2;
+
+ if (generate_breakvec1)
+ {
+ breakvec1.resize(0); // we might be reusing this vector from a previous child
+
+ chromosome->DrawBreakpoints(pedigree_parent1, strand1, strand2, /* p_num_breakpoints */ -1, breakvec1, &heteroduplex1, "addRecombinant()");
+ }
+ else
+ {
+ std::sort(breakvec1.begin(), breakvec1.end());
+ breakvec1.erase(unique(breakvec1.begin(), breakvec1.end()), breakvec1.end());
+ }
+
+ if (generate_breakvec2)
+ {
+ breakvec2.resize(0); // we might be reusing this vector from a previous child
+
+ chromosome->DrawBreakpoints(pedigree_parent2, strand3, strand4, /* p_num_breakpoints */ -1, breakvec2, &heteroduplex2, "addRecombinant()");
+ }
+ else
+ {
+ std::sort(breakvec2.begin(), breakvec2.end());
+ breakvec2.erase(unique(breakvec2.begin(), breakvec2.end()), breakvec2.end());
+ }
+
+ if (breakvec1.size() && (breakvec1.back() > chromosome->last_position_))
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): breaks1 contained a value (" << breakvec1.back() << ") that lies beyond the end of the chromosome." << EidosTerminate();
+ if (breakvec2.size() && (breakvec2.back() > chromosome->last_position_))
+ EIDOS_TERMINATION << "ERROR (Subpopulation::ExecuteMethod_addRecombinant): breaks2 contained a value (" << breakvec2.back() << ") that lies beyond the end of the chromosome." << EidosTerminate();
+
+ // We used to need to look for a leading 0 in the breaks vectors, and swap the corresponding strands,
+ // but HaplosomeRecombined() now handles that for us since it is shared functionality.
+
// Make the new individual as a candidate
Individual *individual = NewSubpopIndividual(/* index */ -1, child_sex, /* age */ 0, /* fitness */ NAN, mean_parent_age);
Haplosome *haplosome1 = haplosome1_null ? chromosome->NewHaplosome_NULL(individual, 0) : chromosome->NewHaplosome_NONNULL(individual, 0);
Haplosome *haplosome2 = nullptr;
if (make_second_haplosome)
- haplosome2_null ? chromosome->NewHaplosome_NULL(individual, 1) : chromosome->NewHaplosome_NONNULL(individual, 1);
+ haplosome2 = haplosome2_null ? chromosome->NewHaplosome_NULL(individual, 1) : chromosome->NewHaplosome_NONNULL(individual, 1);
if (pedigrees_enabled)
{
@@ -7955,6 +8021,9 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
else
{
(population_.*(population_.HaplosomeRecombined_TEMPLATED))(*chromosome, *haplosome1, strand1, strand2, breakvec1, mutation_callbacks);
+
+ if (heteroduplex1.size() > 0)
+ population_.DoHeteroduplexRepair(heteroduplex1, breakvec1, strand1, strand2, haplosome1);
}
}
else
@@ -7999,6 +8068,9 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
else
{
(population_.*(population_.HaplosomeRecombined_TEMPLATED))(*chromosome, *haplosome2, strand3, strand4, breakvec2, mutation_callbacks);
+
+ if (heteroduplex2.size() > 0)
+ population_.DoHeteroduplexRepair(heteroduplex2, breakvec2, strand3, strand4, haplosome2);
}
}
else
@@ -8015,9 +8087,9 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
}
}
}
- else
+ else if (make_second_haplosome)
{
- // both strands are NULL, so we make a null haplosome; we do nothing but record it
+ // both strands are NULL and this is a diploid chromosome, so we record a null haplosome
if (species_.RecordingTreeSequence())
species_.RecordNewHaplosome(nullptr, haplosome2, nullptr, nullptr);
@@ -8030,8 +8102,6 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
chromosome->StopMutationRunExperimentClock("addMultiRecombinant()");
// Run the candidate past modifyChild() callbacks; the target subpop's registered callbacks are used
- bool proposed_child_accepted = true;
-
if (registered_modify_child_callbacks_.size())
{
// BCH 4/5/2022: When removing excess pseudo-parameters from callbacks, we lost a bit of
@@ -8040,7 +8110,10 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
// documented, and I doubt anybody was using it, and they can do the same without the
// modifyChild() callback, so I'm not viewing this loss of functionality as an obstacle
// to making this change.
- proposed_child_accepted = population_.ApplyModifyChildCallbacks(individual, /* p_parent1 */ nullptr, /* p_parent2 */ nullptr, /* p_is_selfing */ false, /* p_is_cloning */ false, /* p_target_subpop */ this, /* p_source_subpop */ nullptr, registered_modify_child_callbacks_);
+ // BCH 1/10/2025: Note that we pass nullptr for the parents, and false for p_is_selfing and
+ // p_is_cloning. We could use pedigree_parent1 and pedigree_parent2, and draw inferences from
+ // them about selfing/cloning, but that would not always be correct; it would just be a guess.
+ bool proposed_child_accepted = population_.ApplyModifyChildCallbacks(individual, /* p_parent1 */ nullptr, /* p_parent2 */ nullptr, /* p_is_selfing */ false, /* p_is_cloning */ false, /* p_target_subpop */ this, /* p_source_subpop */ nullptr, registered_modify_child_callbacks_);
if (!proposed_child_accepted)
{
@@ -8064,35 +8137,37 @@ EidosValue_SP Subpopulation::ExecuteMethod_addRecombinant(EidosGlobalStringID p_
}
}
- if (proposed_child_accepted)
+ if (individual)
{
nonWF_offspring_individuals_.emplace_back(individual);
result->push_object_element_NORR(individual);
#if defined(SLIMGUI)
- gui_offspring_crossed_++;
-
- // This offspring came from parents in various subpops but ended up here, so it is, in effect,
- // a migrant; we tally things, for SLiMgui display purposes, as if it were generated in the
- // parental subpops and then moved. This is gross, but runs only in SLiMgui, so whatever :->
- // We do not set its migrant_ flag, though; that flag is only for takeMigrants() in nonWF.
- // Note that we use pedigree_parent1 and pedigree_parent2 for this, rather than the parents
- // of the strands; in the addMultiRecombinant() case there are potentially many strands with
- // many different parents. It is just too complex to try to keep track of just for SLiMgui.
- // BCH 1/7/2025: It used to be that addRecombinant() did this by strand, but I've dumbed it
- // down to match addMultiRecombinant(); nobody will ever notice. Doing it by strand was
- // approximate too, and maybe even less good in some scenarios.
- if (pedigree_parent1)
- {
- pedigree_parent1->subpopulation_->gui_premigration_size_ += 0.5;
- if (pedigree_parent1->subpopulation_ != this)
- gui_migrants_[pedigree_parent1->subpopulation_->subpopulation_id_] += 0.5;
- }
- if (pedigree_parent2)
- {
- pedigree_parent2->subpopulation_->gui_premigration_size_ += 0.5;
- if (pedigree_parent2->subpopulation_ != this)
- gui_migrants_[pedigree_parent2->subpopulation_->subpopulation_id_] += 0.5;
+ {
+ gui_offspring_crossed_++;
+
+ // This offspring came from parents in various subpops but ended up here, so it is, in effect,
+ // a migrant; we tally things, for SLiMgui display purposes, as if it were generated in the
+ // parental subpops and then moved. This is gross, but runs only in SLiMgui, so whatever :->
+ // We do not set its migrant_ flag, though; that flag is only for takeMigrants() in nonWF.
+ // Note that we use pedigree_parent1 and pedigree_parent2 for this, rather than the parents
+ // of the strands; in the addMultiRecombinant() case there are potentially many strands with
+ // many different parents. It is just too complex to try to keep track of just for SLiMgui.
+ // BCH 1/7/2025: It used to be that addRecombinant() did this by strand, but I've dumbed it
+ // down to match addMultiRecombinant(); nobody will ever notice. Doing it by strand was
+ // approximate too, and maybe even less good in some scenarios.
+ if (pedigree_parent1)
+ {
+ pedigree_parent1->subpopulation_->gui_premigration_size_ += 0.5;
+ if (pedigree_parent1->subpopulation_ != this)
+ gui_migrants_[pedigree_parent1->subpopulation_->subpopulation_id_] += 0.5;
+ }
+ if (pedigree_parent2)
+ {
+ pedigree_parent2->subpopulation_->gui_premigration_size_ += 0.5;
+ if (pedigree_parent2->subpopulation_ != this)
+ gui_migrants_[pedigree_parent2->subpopulation_->subpopulation_id_] += 0.5;
+ }
}
#endif
}