Explain r_break = math.log(2)

[richard-Hobart]

Hi everyone,
Sorry to sound a bit thick on this one, but can anyone perhaps elaborate on why we set the recombination rate at a chromosome break point to equal log(2) ?

the relevant paragraph in
https://tskit.dev/msprime/docs/latest/ancestry.html#example-2-simulating-without-a-pedigree
says

"the probability of co-inheritance across generations of two adjacent base pairs separated by a recombination rate r per base pair is exp(-rt). To make these agree, we set the recombination rate in the base pair segments separating chromosomes to log(2)."

Earlier in the section it says "the probability that two chromosomes are co-inherited across generations is 2^{-t}".

I am just not seeing the connection and why r needs to be bumped up to 0.69315

[petrelharp]

I think the confusion comes from what r is? For msprime, r is not a probability-per-generation, but an (continuous-time) rate. So, the number of crossovers in one unit of time is Poisson with mean r. So, the probability of there being at least one crossover in one unit of time is 1 - exp(-r), which is approximately r - r^2/2. For very small r the difference between "r is the probability of one crossover in one unit of time" and "there are a Poisson(r) number of crossovers per unit time" is neglegible, but when we want the probability of crossover to be 1/2, then, well, it turns out we want r = log(2).

There's a lengthier discussion of related things (although the detalis are different) in Section 14.13 of the SLiM manual.

[richard-Hobart]

Thanks Peter,
the section in the SLiM manual was very helpful. I can now see why we need to transform the requested p value into a repararmeterised λ value using the formula λ = -log(1-p).

I was wondering if you could help me with the following bit of text

"This means that we can add up the λ values for each region along a whole chromosome, even with a complex
recombination map, and draw a preliminary number of crossovers from a Poisson distribution with
the total λ, and then after we select locations with the usual weighted uniform draws, and sort and
unique them, the probability of crossover at every specific site will be as the user requested"

Suppose I had the following rateMap

┌────────────────────────────────────────────────────┐
│left │right │ mid│ span│ rate│
├────────────────────────────────────────────────────┤
│0 │9999999 │ 4999999.5│ 9999999│ 2e-08│
│9999999 │10000000 │ 9999999.5│ 1│ 0.69│
│10000000 │19999999 │ 14999999.5│ 9999999│ 2e-08│
└────────────────────────────────────────────────────┘

does this mean the total λ value = 9999999 * 2e-08 + 1 * 0.69 + 9999999 * 2e-08 = 1.09
I am not sure how SLiM then selects locations with weighted uniform draws. I could not find any specifics in the manual. Perhaps I need to contact the authors

Richard