Is it possible to control the number of SegSites created by MaCS?

[alanaselli]

If I understand correctly, we must provide segSites = NULL to RunMacs or RunMacs2 in order to retrieve all the segregating sites generated. However, in my experience the number of sites generated varies each time I run the command. Furthermore, I would be interested in generating a higher number of sites, if possible. Thus the question: can the user specify the number of sites generated?

Thank you!

[gregorgorjanc]

Yes, the number of generated sites varies between simulation runs - it’s an outcome of a stochastic process.

If you want more segregating sites, then the implemented demographic model for the species is not matching your expectation. Have a look at runMacs2() if you want to define demographic model yourself.

[alanaselli]

Thank you very much for the quick response, that solves my problem.

[gaynorr]

There are several parameters that will affect the number of sites that are simulated. The most important are the number of individuals, the effective population size (Ne), and the historic bottleneck(s) that are represented using the histNe and histGen arguments. The number of base pairs, genetic length, and mutation rate also matter, but we should have reasonable estimates for these values.

Here is a little illustration of how bottlenecks effect the number of sites that get simulated.

library(AlphaSimR)

# No historic bottleneck, current population modeled as having Ne 100 
# Runs very quickly
A = runMacs2(nInd = 100,
             nChr = 1, 
             Ne = 100, # Current effective population size
             histNe = NULL,
             histGen = NULL)
A@nLoci # ~6,000 loci are simulated

# Add a bottleneck 100 generations ago. The populations goes from Ne 1000 to 100.
# Takes a while longer to run
B = runMacs2(nInd = 100,
             nChr = 1, 
             Ne = 100, # Current effective population size
             histNe = 1000, # Historic Ne of 1000
             histGen = 100) # Number of generations ago for historic Ne
B@nLoci # ~22,000 loci are simulated

# The bottleneck is changed to reflect Ne going from 10,000 to 100.
# Takes even longer to run
C = runMacs2(nInd = 100,
             nChr = 1, 
             Ne = 100, # Current effective population size
             histNe = 10000, # Historic Ne of 1000
             histGen = 100) # Number of generations ago for historic Ne
C@nLoci # ~190,000 loci are simulated

[alanaselli]

Thank you @gaynorr, this gives further understanding of how the parameters affect the simulation.
In my case, I'm simulating a Nelore (cattle) population, so it makes sense to use the argument species. I was experimenting simulating the entire history on AlphaSimR (based on a study that used another Software), and trying to make it resemble more Nelore, which is a bit different than taurine breeds.

[gregorgorjanc]

@alanaselli the implemented cattle demographic model is based on Holstein. You might want to evaluate if you need a different demographic model.