Using BDMM with MLVA 15 data (BEASTvntr)

[Sheldon Waugh]

Good evening, 

I'm trying to run a model using the multi-type Birth-Death model combined with the a Sainudiin substitution model based on VNTR data. I wouldn't think there would be any serious problems, however I encountered a fatal error when running this on Beast 2.4.1.

java.lang.ArrayIndexOutOfBoundsException: 0

        at org.ejml.data.DenseMatrix64F.<init>(DenseMatrix64F.java:97)

        at beast.evolution.substitutionmodel.EJMLEigenSystem.decomposeMatrix(Unknown Source)

        at beast.evolution.substitutionmodel.Sainudiin.getEigenDecomposition(Unknown Source)

        at beast.evolution.likelihood.BeagleTreeLikelihood.setUpSubstModel(Unknown Source)

        at beast.evolution.likelihood.BeagleTreeLikelihood.initialize(Unknown Source)

        at beast.evolution.likelihood.BeagleTreeLikelihood.initAndValidate(Unknown Source)

        at beast.core.BEASTInterface.initByName(Unknown Source)

        at beast.evolution.likelihood.TreeLikelihood.initAndValidate(Unknown Source)

        at beast.util.XMLParser.initBEASTObjects(Unknown Source)

        at beast.util.XMLParser.parse(Unknown Source)

        at beast.util.XMLParser.parseFile(Unknown Source)

        at beast.app.BeastMCMC.parseArgs(Unknown Source)

        at beast.app.beastapp.BeastMain.main(Unknown Source)

        at beast.app.beastapp.BeastLauncher.main(Unknown Source)

TreeLikelihood(treeLikelihood.KZ-Brucella-VNTR-15-final_w_spatial_081417) uses BeerLikelihoodCore

  Alignment(KZ-Brucella-VNTR-15-final_w_spatial_081417): [taxa, patterns, sites] = [264, 15, 15]

Type trait with the following types detected:

Animal (0)

Human (1)

===============================================================================

Citations for this model:

Bouckaert RR, Heled J, Kuehnert D, Vaughan TG, Wu C-H, Xie D, Suchard MA,

  Rambaut A, Drummond AJ (2014) BEAST 2: A software platform for Bayesian

  evolutionary analysis. PLoS Computational Biology 10(4): e1003537

Kuehnert D, Stadler T, Vaughan TG, Drummond AJ. 2016. Phylodynamics with migration:

        A computational framework to quantify population structure from genomic data.

        Mol Biol Evol. 33(8):2102–2116.

Timothy G. Vaughan, Denise Kuhnert, Alex Popinga, David Welch and

Alexei J. Drummond, 'Efficient Bayesian inference under the

structured coalescent', Bioinformatics 30:2272, 2014.

Raazesh Sainudiin et al. (2004) Microsatellite Mutation Models.

  Genetics 168:383-395

Chieh-Hsi Wu and  Alexei J. Drummond. (2011) Joint Inference of

  Microsatellite Mutation Models, Population History and Genealogies

  Using Transdimensional Markov Chain Monte Carlo.

  Genetics 188:151-164

===============================================================================

Start likelihood: -Infinity after 10 initialisation attempts

P(posterior) = -Infinity (was -Infinity)

        P(prior) = -Infinity (was -Infinity)

                P(birthDeathMigration.t:KZ-Brucella-VNTR-15-final_w_spatial_081417) = -Infinity (was -Infinity)

                P(RPrior.t:KZ-Brucella-VNTR-15-final_w_spatial_081417) = -3.7046244414474376 (was -3.7046244414474376)

                P(a1Prior.s:KZ-Brucella-VNTR-15-final_w_spatial_081417) = 0.0 (was 0.0)

                P(b0Prior.s:KZ-Brucella-VNTR-15-final_w_spatial_081417) = 0.0 (was 0.0)

                P(b1Prior.s:KZ-Brucella-VNTR-15-final_w_spatial_081417) = 0.0 (was 0.0)

                P(becomeUninfectiousRatePrior.t:KZ-Brucella-VNTR-15-final_w_spatial_081417) = -1.8378770664093453 (was -1.8378770664093453)

                P(ClockPrior.c:KZ-Brucella-VNTR-15-final_w_spatial_081417) = 0.0 (was 0.0)

                P(gPrior.s:KZ-Brucella-VNTR-15-final_w_spatial_081417) = 0.0 (was 0.0)

                P(rateMatrixPrior.t:KZ-Brucella-VNTR-15-final_w_spatial_081417) = -0.20000000000000012 (was -0.20000000000000012)

                P(samplingProportionPrior.t:KZ-Brucella-VNTR-15-final_w_spatial_081417) = 1.7763568394002497E-15 (was 1.7763568394002497E-15)

        P(likelihood) = NaN (was NaN)  **

                P(treeLikelihood.KZ-Brucella-VNTR-15-final_w_spatial_081417) = NaN (was NaN)  **

Fatal exception: Could not find a proper state to initialise. Perhaps try another seed.

java.lang.RuntimeException: Could not find a proper state to initialise. Perhaps try another seed.

        at beast.core.MCMC.run(Unknown Source)

        at beast.app.BeastMCMC.run(Unknown Source)

        at beast.app.beastapp.BeastMain.<init>(Unknown Source)

        at beast.app.beastapp.BeastMain.main(Unknown Source)

        at beast.app.beastapp.BeastLauncher.main(Unknown Source)

Fatal exception: Could not find a proper state to initialise. Perhaps try another seed.

java.lang.RuntimeException: An error was encounted. Terminating BEAST

        at beast.app.util.ErrorLogHandler.publish(Unknown Source)

        at java.util.logging.Logger.log(Logger.java:738)

        at java.util.logging.Logger.doLog(Logger.java:765)

        at java.util.logging.Logger.log(Logger.java:788)

        at java.util.logging.Logger.severe(Logger.java:1463)

        at beast.app.beastapp.BeastMain.<init>(Unknown Source)

        at beast.app.beastapp.BeastMain.main(Unknown Source)

        at beast.app.beastapp.BeastLauncher.main(Unknown Source)

I am attaching the XML file to this post.

I appreciate the help in advance!!

-Sheldon

[Tim Vaughan]

Hi Sheldon, if at all possible, I'd try updating the version of BEAST 2 you're using - 2.4.7 is the latest release.  If the problem remains, it looks as though this is to do with the substitution model you've selected, as it's the tree likelihood that's producing the error.

Tim

--
You received this message because you are subscribed to the Google Groups "beast-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to beast-users...@googlegroups.com.
To post to this group, send email to beast...@googlegroups.com.
Visit this group at https://groups.google.com/group/beast-users.
For more options, visit https://groups.google.com/d/optout.

[Sheldon Waugh]

Thanks for the advice Tim, 

No dice on 2.4.7.

I guess I can't combine this substitution model with the BDMM package.

I don't know why, though.

I appreciate it. 

Thanks,

Sheldon

[Tim Vaughan]

Hi again Sheldon,

I've had another look at this and I don't think BDMM is the issue here. There _is_ a problem with the BDMM side of your analysis setup.  (This is causing the error you see at the end of the output, where beast is complaining that it can't initialize the chain.) It seems there are some samples in your dataset that occur during periods where you've set the sampling proportion to 0:  the samplingProportion.t parameter in your XML has zeros in it, and replacing these with non-zero values allows the BDMM tree prior to produce a non-zero density on the starting state.

However, this doesn't fix the other problem, which is with the tree likelihood calculation under the VNTR model. The problem occurs even when I comment out the BDMM portion of your model, so I doubt it's a problem with the interaction between the packages.  I've never used this model before though, so hopefully someone else will have some suggestions on how to solve this.

Tim

[Denise]

Hi Sheldon,

I agree with Tim. You might have to post a new question where the subject line states clearly that you have a problem with this particular substitution model, so people who've used it can help you out. 

For BDMM you might want to use multi-rho-sampling instead of serial sampling (see section 2.2.3 of the bdsky tutorial: https://github.com/BEAST2-Dev/bdsky/blob/master/doc/bdsky_tutorial_2.1.1.pdf) because you have many samples at once in each of the sampling years. 

(As a quick fix you can adjust the samplingRateChangeTimes to include all samples like so:

                <parameter id="samplingRateChangeTimes.s:KZ-Brucella-VNTR-15-final_w_spatial_081417" dimension="2" name="samplingRateChangeTimes">0.0 6.01</parameter> )

Denise

[Sheldon Waugh]

Thanks Denise and Tim, 

I will go ahead and create a new question looking at the proper use of my substitution model.

To Denise: I think the part that confused me is when I chose the multi-rho sampling, BEAUTI didn't seem to work correctly, it may have been a mistake on my part, but likely it maybe due to the fact that I have to go into the actual XML and fix it myself.

Additionally, I guess I'll have to get more acclimated with the XML file format!

Thanks again, 

Sheldon Waugh