Estimation from SFS; tpl error extinct deme

[Mitch Lky]

Hello and thank you for any help!

I am trying to perform demographic inference from the SFS and am running into some issues with figuring out the tpl file. I attach a diagram of the models I am trying to define in the tpl and est files.



As you can see it is a 3-4 epoch model of a population decline and subsequent 4-5 way divergence. When I run fsc28 with the attached tpl and est files I get the following ;
> fsc28 -t ${name}.tpl -n100000 --msfs -e ${name}.est -M -L40 -c0 -B10

> Error: deme 0 going extinct while still having remaining lineages

Looking over somethings online I tried redefining things where instead of resizing as a complex parameter I used absoluteResize. I also tried recoding the coalescence of N0 into NBOT as TDIV 0 0 1, because one example online did so but migrating from deme 0 to deme 0 doesn't make any sense. Neither of these tests changed the error.

I would love if someone could look over my tpl and est files and help me see what I am missing.

Thank you for your help!
Mitch Lokey
Cornell University

[Mitch Lky]

In case it is helpful to not download the files I post the 3 epoch files here.

era3-mig0.tpl
//Parameters for the coalescence simulation program : fastsimcoal.exe
5 samples to simulate :
//Population effective sizes (number of genes haploid)
NABS
NONF
NJDSP
NPLE
NSPSP
//Samples sizes, samples age, inbreeding ONF=0.12 , ABS=0.155, JDSP=0.17, SPSP=0.15
48
48
48
48
48
//Growth rates : negative growth implies population expansion
0
0
0
0
0
//Number of migration matrices : 0 implies no migration between demes
0
//historical event: time, source, sink, migrants, new deme size, new growth rate, migration matrix index
6 historical event
TDIV 1 0 1 0 0 0
TDIV 2 0 1 0 0 0
TDIV 3 0 1 0 0 0
TDIV 4 0 1 0 0 0
TDIV 0 0 0 BOTSIZE 0 0
TBOT 0 0 0 ANCSIZE 0 0
// Also tried with; TDIV 0 0 1 NBOT 0 0 absoluteResize
//Number of independent loci [chromosome]
1 0
//Per chromosome: Number of linkage blocks
1
//per Block: data type, num loci, rec. rate and mut rate + optional parameters
FREQ 1000000 1.5e-8 3e-9 0

era3-mig0.est
// Priors and rules file
// *********************

[PARAMETERS]
//#isInt? #name #dist.#min #max
//all Ns are in number of haploid individuals
1 NANC      unif 1   100000   output
1 NBOT      unif 1   100000   output
1 NONF      unif 1   100000   output
1 NABS      unif 1   100000   output
1 NPLE      unif 1   100000   output
1 NSPSP     unif 1   100000   output
1 NJDSP     unif 1   100000   output
1 TDIV      unif 1   100000   output
1 TPLUSDIV  unif 1   100000   output

[COMPLEX PARAMETERS]
0 ANCSIZE = NANC/NBOT       output
0 BOTSIZE = NBOT/NONF       output
1 TBOT = TDIV+TPLUSDIV      output

Thanks again!

[Mitch Lky]

Hello everyone,

I solved the error (thanks Elissa). It turned out I wasn't thinking about migration and new deme size for the populations divergence correctly.

It ran by editing the above code to;
TDIV 1 0 1 1 0 0
TDIV 2 0 1 1 0 0
TDIV 3 0 1 1 0 0
TDIV 4 0 1 1 0 0

TDIV 0 0 0 BOTSIZE 0 0

However, the estimates I got are pretty wild and far from our estimates from SMC++ and dadi.

Can anyone confirm that my tpl is setting up the model the way I intend in the attached figure?

Thanks!
Mitch