@jeromekelleher
Please check this: the test test_gaps in test models fails for the smc_prime (now mapped to smc(k=1)). This tests overlapping segments in a non-discrete genome. I think smc(k=1) is not equavelent to the smc_prim in the case of a non-discrete genome. What do you think?

Please check this: the test test_gaps in test models fails for the smc_prime (now mapped to smc(k=1)). This tests overlapping segments in a non-discrete genome. I think smc(k=1) is not equavelent to the smc_prim in the case of a non-discrete genome. What do you think?

Ooh - that's a nasty one! Can you confirm by making k=1.00001 and rerunning?

ts = msprime.simulate(
                sample_size=sample_size,
                recombination_rate=recombination_rate,
                random_seed=random_seed,
                model=msprime.SmcKApproxCoalescent(hull_offset=1.00001),
            )
            edgesets = sorted(ts.edgesets(), key=lambda e: (e.parent, e.left))
            num_found = 0
            for j in range(1, len(edgesets)):
                r = edgesets[j - 1]
                s = edgesets[j]
                if r.right != s.left and r.parent == s.parent:
                    num_found += 1
            assert num_found == 0 #fails

hull_offset=1.00001 gives similar results to 1. both fail the test

Hmm, good catch. We'll have to think about that.

Ah right, the problem is k is too big. Here, the genome is one unit long, so k=1 means we allow all possible overlaps. If we put k=0.00001 in this example we get the right answer.

Now, I'm not entirely sure what the right value of the smc_prime is in a continuous genome model. I guess it's a very small number? It's an annoying one.

In a math term, it should be when k -> 0. I am not sure if we set it to the smallest float that could be stored, it would be a sustainable universal solution. Also, we would have to set the k value differently for discrete genomes (k=1).

Let's try a small value like 1e-14 or something so as the default, and see if that does the right thing for the discrete genome case. Having a different mapping for discrete genome and not would be quite tricky, as we'd have to pass the discrete_genome parameter around a lot, and it would be confusing to document anyway.

Using very small values like flt_min or whatever could have unexpected consequences, so I think a reasonable like 1e-14 is less likely overall to bite us.

It does not break unit tests, and I am checking the stats verifications. Verifications does not look ok though. Check below the number of trees between msprime.SmcKApproxCoalescent(hull_offset=1e-14) and msprime.SmcPrimeApproxCoalescent()

for: the test_smck_vs_smckapprox_recomb_discrete_hotspots

recombination_map = msprime.RateMap(
           position=[0, 100, 500, 900, 1200, 1500, 2000],
           rate=[0.00001, 0, 0.0002, 0.00005, 0, 0.001],
       )
       self._run(
           samples=10,
           population_size=1000,
           recombination_rate=recombination_map,
           num_replicates=300,
           discrete_genome=True,
       )

What about we restrict the smc_prime model to discrete genome from the interface?

What about we restrict the smc_prime model to discrete genome from the interface?

I considered this, but it's not trivial because of the complex structure of the model argument and also it will likely break code.

Hmm, this is annoying.

OK, let's drop the remapping of the strings for now, and keep the reparameterisation of smck to require the parameter. At least this excecise has shown us that a default of 1 was definitely not the right thing to do...

It may be simpler to just do this in a different PR or update this one, it's up to you.

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 91.22%. Comparing base (1526582) to head (902179d).
⚠️ Report is 2 commits behind head on main.

@@           Coverage Diff           @@
##             main    #2445   +/-   ##
=======================================
  Coverage   91.22%   91.22%           
=======================================
  Files          20       20           
  Lines       11847    11847           
  Branches     2300     2300           
=======================================
  Hits        10807    10807           
  Misses        569      569           
  Partials      471      471           

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

check now. Also check the recommendation in the documentation: for smc_prime the documentation recommends smck(1) for discrete and smck(1e-14) for continuous genomes

Docstring is malformed, that's why the docs build is failing

yes, should be fixed now

@jeromekelleher Can be merged