NAN confidence intervals FIM

[Zoe Meziere]

Hi Ryan and dadi users,

I am trying to get confidence intervals using Godambe.FIM_uncert() but I am running into this error at every line:

WARNING:Inference:Model is masked in some entries where data is not.

which result in this:

Estimate parameter standard deviations from FIM with eps=0.1: [nan nan nan nan nan nan]

Do you know what could be causing the issue?

Thank you in advance for your help,

Best,

Zoe

[Ryan Gutenkunst]

Hello Zoe,

That Warning suggests that dadi is having trouble calculating the model frequency spectrum for at least some set of parameters that is used in the FIM calculation. This might result from parameters that cross zero when calculating the numerical derivatives. I would suggest trying a smaller eps setting and seeing whether the error persists.

Best,
Ryan

> --
> You received this message because you are subscribed to the Google Groups "dadi-user" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to dadi-user+...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/dadi-user/3eca7b38-89a6-46ce-9400-58bd3d2d5c03n%40googlegroups.com.

[Zoe Meziere]

Thanks for your quick reply Ryan!

I have tried smaller EPS settings and keep having the same issue and if I go too low, I get the following error message:

WARNING:Numerics:Extrapolation may have failed. Check resulting frequency spectrum for unexpected results.
Traceback (most recent call last):
  File "FIM_confindence_intervals.py", line 88, in <module>
    main(snps, model, sims, eps, opt, PTS)
  File "FIM_confindence_intervals.py", line 52, in main
    ll_model = Inference.ll_multinom(sim_model, fs)
  File "/home/uqzmezie/.conda/envs/dadi/lib/python3.7/site-packages/dadi/Inference.py", line 563, in ll_multinom
    ll_arr = ll_multinom_per_bin(model, data)
  File "/home/uqzmezie/.conda/envs/dadi/lib/python3.7/site-packages/dadi/Inference.py", line 547, in ll_multinom_per_bin
    return ll_per_bin(theta_opt*model, data)
  File "/home/uqzmezie/.conda/envs/dadi/lib/python3.7/site-packages/dadi/Inference.py", line 484, in ll_per_bin
    if hasattr(data, 'folded') and data.folded and not model.folded:
AttributeError: 'MaskedArray' object has no attribute 'folded'

Do you think it could have to do with my bootstraps? Should I re-run the boostraping with smaller EPS settings?

Many thanks,

Zoe

[Ryan Gutenkunst]

Hello Zoe,

The bootstraps aren’t affect by eps, so it would surprise me if that was the issue. And the FIM doesn’t require bootstraps. Are you running this through dadi-cli, or the Python interface?

Best,
Ryan

> To view this discussion on the web visit https://groups.google.com/d/msgid/dadi-user/99c4c58b-b141-4328-8225-7c55d3aaf237n%40googlegroups.com.

[Zoe Meziere]

Hi Ryan, and thank you for your help! I am using the python interface.

I am been trying to get to the root of the issue and plotted the data and modelled SF spectra, and got these weird plots (see Plot1 attached). We think this might be because the model is unrealistic. Indeed, I have been having issues with convergence and parameters hitting the set boundaries. 

These populations are very little differentiated - Fst = 0.01 and panmictic looking on PCA. I used the scramble_pop_ids() function, and from the plots obtained (see attached Plot2, Plot3, Plot4), I understand they are indeed little differentiated. Is my interpretation correct?

To give you some context, I am using a simple asymmetric migration model to estimate migration rates and gene flow. I am not interested in comparing different demographic models. Do you think migration rates and/or population sizes might be too high here? Any any idea how to get around this issue?

Thank you again!!! 

[Ryan Gutenkunst]

Hello Zoe,

Yes, these populations look almost panmictic. That would lead to issues with model fitting, since if migration rates get very high, the numerical solution breaks down. In general, we never set migration rates about 10 or 20, since that leads to something like panmixia anyways. The easiest way to address numerical issues is to increase your calculation grid points. That makes it more accurate at the cost of slower calculation.

Best,
Ryan

> To view this discussion on the web visit https://groups.google.com/d/msgid/dadi-user/18d33af1-6b4b-4389-940f-db876cae17c2n%40googlegroups.com.
> <Plot2_scrambledSFS.pdf><Plot3_dataSFS.pdf><Plot1_modelledVSdata.pdf><Plot4_residualsSFS.pdf>

[Zoe Meziere]

Hi Ryan,

Thank you for your help!

For a population pair with projection [70,20], what do you think would be good grid points? I have tried PTS=[80,90,100]. Would you increase it to [100,110,120] in this case?

And do you think these parameter bounds are sensible?

migration rates: 0.001 - 20

population sizes: 0.001 - 150

divergence time: 0.001 - 5

Thanks again,

Zoe

[Ryan Gutenkunst]

Hello Zoe,

You can increase the grid size to as much as you want. Trying [100,110,120] would be reasonable, but if you had the compute power you could go to [200,210,220] easily.

Those boundaries are pretty typical, although you can easily take the lower boundary of the divergence time and migration rates to 0. (Those aren’t numerically problematic.)

Best,
Ryan

> To view this discussion on the web visit https://groups.google.com/d/msgid/dadi-user/56ac6136-04b9-48cc-a235-4916a6249a16n%40googlegroups.com.

[Zoe Meziere]

Great, thank you very much for you help Ryan!

Best, 

Zoe