read-in custom ratemap

[archgen]

Hi, just checking if it is possible to read in my own ratemap.txt file in the ratemap format of (left \t right \t mid \t span \t rate)? I see there is the read_hapmap() but wondering if msprime can read-in a custom ratemap.
Thanks very much!!

[grahamgower]

Hey @archgen! 👋

This should be possible, assuming your rate column uses the units cM/Mb (if not, you'll have to convert), and assuming your map intervals are contiguous (otherwise you'll have to figure out what rate to give the "gaps"). Something like the following should work.

chromosome_length = ... # some appropriate number here
ratemap = msprime.RateMap.read_hapmap(
    "ratemap.txt", position_col=0, rate_col=4, sequence_length=chromosome_length
)

If you get errors from this, you might need to manually edit your file (e.g. add entries to ensure the telomeres have rate zero). See the API docs for info on the meaning of the columns.
https://tskit.dev/msprime/docs/stable/api.html#msprime.RateMap

And then see the discussion of rate maps to see how to inspect the data structure once you've read it in, to double check that you got the ratemap you expect.
https://tskit.dev/msprime/docs/stable/rate_maps.html