Stdpopsim: add first background selection (BGS) model

[grahamgower]

Now that #560 has been merged, we should choose a model from #391 and implement it (using #587, once its done). This first model will help us think about what other infrastructure is needed, like what form the BGS catalog will take.

[andrewkern]

@mufernando ping me if you want to move the internals of the annotations over to pyranges. should be very easy if we need this.

[mufernando]

Thanks! I think we can leave annotations as it is right now, but use pyranges to do operations with the new GenomicElements class.

I'll get started on implementing this today/wknd.

[mufernando]

One thing we haven't considered yet is that those DFEs are tied to specific demographic models (which may not have been implemented in stdpopsim yet).

Do we want to make DFEs property of a demographic model or not?

[apragsdale]

Do we want to make DFEs property of a demographic model or not?

I think that would end up being too restrictive. Most inferred DFEs use data from a single population and thus would be constrained to single-population settings. You wouldn't be able to run, e.g., the out of Africa demography with a DFE in that case.

[grahamgower]

I agree with Aaron here. And I'm hopeful that the ability to simulate BGS models easily will accelerate the development of more generally applicable DFEs. (Or conversely, show that we really should be more restrictive)

[jeromekelleher]

Also +1 for decoupling demography from other evolutionary processes.

[mufernando]

I agree, but we need to think about how to decouple the DFE from demography.

Nes is inferred in these papers so:

• Should we get s and add that to the catalog? If so, it's not necessarily clear how to go from Nes to s because Ne changes with time and there are multiple pops.
• Or should we add Nes to the catalog, figure out a way to calculate Ne from a given demographic model and then use that to establish the s?

[grahamgower]

I think we have related problems with our use of Ne and mutation rate. A paper infers Ne with specific assumptions on mutation rate and generation time, and then along we come and simulate that Ne with a different mu. See #557. I think there was some more discussion of this elsewhere, but I can't seem to find it right now. I think the consensus was "that's a hard problem to deal with in a pragmatic way, we look forward to seeing the solution!"

[mufernando]

yeah, definitely something we should explore a bit with analyses for the next paper.

I think one practical way forward would be to save DFEs to the catalog using the Ne*s parametrization, given that is what methods actually infer. But, when choosing to simulate from a particular DFE, we can get at s by using @petrelharp coalescence rate calculation (which was added to msprime).

[grahamgower]

I reckon we really do want to specify s for simulations, not Ne*s, because the effect of purifying selection should be more intense with larger populations, as weakly deleterious alleles will be removed more efficiently in bigger populations. If we adjusted s based on the current Ne, then the effect of BGS would be identical for all values of Ne.

So this leaves open the question of what value of s to use in the event that Ne*s has been reported, rather than s itself. This is outside my expertise, so it would be nice if others weighed in... Maybe use the harmonic mean of Ne(t), where Ne(t) is the model used for inference of Ne*s (probably constant anyway)?

[apragsdale]

I agree with @grahamgower that we want to have the distributions of s instead of 2*Ne*s. I think most DFE papers that I've read have converted their inferred distributions of gamma and report distributions of s. But in the case that they haven't, if it's inferred in a dadi framework, then typically the Ne refers to the ancestral/reference Ne.

And if no Ne is provided or estimated in some study study (i.e. only have demographic models with relative sizes), we could either use our default Ne for the species that we've implemented, or probably better to just not use that study altogether.

[mufernando]

Aaron, can you point me to the papers you are thinking about?

I was looking at the papers compiled in this issue here: #391. Sorry I forgot to mention this earlier.

[mufernando]

could use the input of @dschride @petrelharp

[apragsdale]

Aaron, can you point me to the papers you are thinking about?

I was looking at the papers compiled in this issue here: #391. Sorry I forgot to mention this earlier.

Ok thanks. Yes, some of those have converted to s, though it looks like some of them haven't. If the papers that have not converted don't report Ne (which would usually require some assumption about the per-base mutation rate), we'd probably have to do that ourselves. In which case, maybe we just take the default rates defined for our species?

[petrelharp]

I agree - deal with s, not Ne * s; if papers do not convert then do something reasonable. We maybe should include a scaling parameter, so people could multiply all the ss by something when using it? As for what "something reasonable" is: one option is "not use that paper". Or... there's not an obvious "right" choice of Ne to use in the case of non-constant population size (as far as I can tell?), so "default Ne for the species", which is roughly mean TMRCA, is as good as anything, I think.