Fixing a starting tree in SNAPP

535 views
Skip to first unread message

bander...@gmail.com

unread,
Jan 25, 2016, 4:54:20 AM1/25/16
to beast-users
Hi BEAST users,

I've been trying to get SNAPP to run with two modifications to the xml created by BEAUTi. The first modification is an ingroup prior (making most of my taxa monophyletic to ensure the outgroup is an outgroup). When I try to run this with one of my data sets, it fails to initialise:

Start likelihood: -Infinity after 10000 initialisation attempts
P(posterior) = -Infinity (was NaN)
P(prior) = -Infinity (was NaN)
P(lambdaPrior.even_nopoly_data) = 4.605170185988091 (was NaN)
P(snapprior.even_nopoly_data) = -909.0275559506407 (was NaN)
P(uPrior.even_nopoly_data) = 0.5276471940468044 (was NaN)
P(vPrior.even_nopoly_data) = -1.1872459785043712 (was NaN)
P(ingroup.prior) = -Infinity (was NaN)
P(ingroup.prior) = -Infinity (was NaN)
P(likelihood) = NaN (was NaN)
P(treeLikelihood.even_nopoly_data) = NaN (was NaN)

I assume this means that the ingroup is not monophyletic in the starting upgma trees?

To try to get around this, I looked at (http://beast2.org/fix-starting-tree/) and attempted to make the adjustments to my xml file. I kept running into errors however, and the only way it seemed to work was when I removed the clock.rate parameter. This seems like an undesirable option, and I don't know how SNAPP uses that parameter. Here's what I did:

ORIGINAL XML

<run id="mcmc" spec="MCMC" chainLength="400000" numInitializationAttempts="10000" storeEvery="400">
    <state id="state" storeEvery="400">
        <stateNode id="Tree.even_nopoly_data" spec="beast.util.ClusterTree" clusterType="upgma" nodetype="snap.NodeData">
            <taxa id="snap.even_nopoly_data" spec="snap.Data" dataType="integerdata">
                <rawdata idref="even_nopoly_data"/>
                <taxonset id="peed" spec="TaxonSet">
                    <taxon id="peed.Stew11" spec="Taxon"/>
                    ...
            </taxa>
            <parameter id="RealParameter.0" lower="0.0" name="clock.rate" upper="0.0">1.0</parameter>
        </stateNode>
        <parameter id="u" lower="0.0" name="stateNode">0.589991473402598</parameter>
        <parameter id="v" lower="0.0" name="stateNode">3.27804096929902</parameter>
        <parameter id="lambda" lower="0.0" name="stateNode">0.01</parameter>
        <parameter id="coalescenceRate" name="stateNode">10.0</parameter>
    </state>

MODIFIED XML

<run id="mcmc" spec="MCMC" chainLength="400000" numInitializationAttempts="10000" storeEvery="400">
    <state id="state" storeEvery="400">
        <stateNode id="Tree.even_nopoly_data" spec="beast.util.TreeParser" IsLabelledNewick="true" adjustTipHeights="false" nodetype="snap.NodeData" newick="((((peed,(shova,(shovb,lan))),bas),(nana,(shovel,war))),con);">
            <taxa id="snap.even_nopoly_data" spec="snap.Data" dataType="integerdata">
                <rawdata idref="even_nopoly_data"/>
                <taxonset id="peed" spec="TaxonSet">
                    <taxon id="peed.Stew11" spec="Taxon"/>
                  ...
            </taxa>
            <parameter id="RealParameter.0" lower="0.0" name="clock.rate" upper="0.0">1.0</parameter>
        </stateNode>
        <parameter id="u" lower="0.0" name="stateNode">0.589991473402598</parameter>
        <parameter id="v" lower="0.0" name="stateNode">3.27804096929902</parameter>
        <parameter id="lambda" lower="0.0" name="stateNode">0.01</parameter>
        <parameter id="coalescenceRate" name="stateNode">10.0</parameter>
    </state>

PRODUCES ERROR:

Error 130 parsing the xml input file

This BEASTInterface (Tree.even_nopoly_data) has no input with name clock.rate. Choose one of these inputs: IsLabelledNewick,taxa,newick,offset,threshold,singlechild,adjustTipHeights,scale,initial,trait,taxonset,nodetype,estimate

Error detected about here:
  <beast>
      <run id='mcmc' spec='MCMC'>
          <state id='state'>
              <stateNode id='Tree.even_nopoly_data' spec='beast.util.TreeParser'>


If I remove the clock.rate parameter line, the run starts without a problem.

Does anyone have a suggestion for properly setting the starting tree without removing the clock.rate parameter? Comments?

Remco Bouckaert

unread,
Jan 25, 2016, 2:25:13 PM1/25/16
to beast...@googlegroups.com
By default, SNAPP assumes a fixed clock rate of 1, so the height of the tree is in units of substitutions per site. You should not need to specify a clock rate anywhere, and the version that ran is probably just fine. 

However, if you add timing information for instance by having an age distribution specified for the ingroup prior, you do need to add the clock rate to the state and clock rate model, and you need to estimate the clock rate and have a prior on that rate as well, so there would be more involved than just changing the starting tree.

Hope this helps,

Remco



--
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.

Benjamin Anderson

unread,
Jan 26, 2016, 1:39:23 AM1/26/16
to beast...@googlegroups.com
Thanks Remco!
I don't have any age information (just topology), so I'll give it a go without the clock.rate setting.
Cheers,
Ben

--
You received this message because you are subscribed to a topic in the Google Groups "beast-users" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/beast-users/7Lhl2mffI1o/unsubscribe.
To unsubscribe from this group and all its topics, send an email to beast-users...@googlegroups.com.

Krzysztof Kozak

unread,
Jan 24, 2019, 9:12:22 AM1/24/19
to beast-users
Dear Remco,

Is there an example of an xml with clock parameters, please? I haven't found it in the Github examples.
I have a SNAP xml that works perfectly fine: however, adding clock parameters leads to errors that surprise me, e.g. "This BEASTInterface (Tree.erato25) has no input with name stateNode." - even though the <state> block has not been changed.

Best,
Reply all
Reply to author
Forward
0 new messages