why was my time result less than the calibration time use mcmctree?

[cqm19...@126.com]

Hello everyone,

      I'm caculating divergence time using mcmctree block of the PAML4.9. However my estimated time was less than the calibration time.

      the output files were attached.

      The input tree was as follows: 

14    1

((((((((Capra_hircus,Capra_aegagrus)'>0.20<0.28',Capra_ibex),Pseudois_nayaur),Ammotragus_lervia),((Ovis_canadensis,Ovis_dalli),((Ovis_aries_musimon,Ovis_aries),Ovis_ammon))),(Budorcas_taxicolor,Naemorhedus_goral)),Pantholops_hodgsonii)'>8.9<13.0',Bos_taurus);

//end of file

     

      The output tree was as follows:

#NEXUS

BEGIN TREES;

        UTREE 1 = ((((((((Capra_hircus: 0.267024, Capra_aegagrus: 0.267024) [&95%={0.225842, 0.305493}]: 0.022642, Capra_ibex: 0.289667) [&95%={0.228111, 0.597342}]: 0.763496, Pseudois_nayaur: 1.053163) [&95%={0.671254, 1.59468}]: 0.370653, Ammotragus_lervia: 1.423816) [&95%={0.965885, 2.13695}]: 0.469252, ((Ovis_canadensis: 0.390740, Ovis_dalli: 0.390740) [&95%={0.0751367, 0.957385}]: 0.570426, ((Ovis_aries_musimon: 0.467613, Ovis_aries: 0.467613) [&95%={0.114136, 1.23938}]: 0.260268, Ovis_ammon: 0.727881) [&95%={0.265375, 1.54298}]: 0.233285) [&95%={0.468601, 1.62205}]: 0.931901) [&95%={1.257, 3.54941}]: 0.575487, (Budorcas_taxicolor: 1.496072, Naemorhedus_goral: 1.496072) [&95%={0.48157, 2.61561}]: 0.972483) [&95%={1.51475, 4.05632}]: 0.498584, Pantholops_hodgsonii: 2.967139) [&95%={1.70009, 5.144}]: 14.703475, Bos_taurus: 17.670614) [&95%={15.9542, 19.9371}];

END;

The mcmctree.ctl file was as follows:

           seed = -1

       seqfile = zz_total.txt

      treefile = zz_total.trees

      mcmcfile = mcmc.txt

       outfile = out.txt

         ndata = 1

       seqtype = 0    * 0: nucleotides; 1:codons; 2:AAs

       usedata = 1    * 0: no data; 1:seq like; 2:normal approximation; 3:out.BV (in.BV)

         clock = 2    * 1: global clock; 2: independent rates; 3: correlated rates

       RootAge = '>16.0<20.0'  * safe constraint on root age, used if no fossil for root.

         model = 4    * 0:JC69, 1:K80, 2:F81, 3:F84, 4:HKY85

         alpha = 0    * alpha for gamma rates at sites

         ncatG = 5    * No. categories in discrete gamma

     cleandata = 0    * remove sites with ambiguity data (1:yes, 0:no)?

       BDparas = 1 1 0.1  * birth, death, sampling

   kappa_gamma = 6 2      * gamma prior for kappa

   alpha_gamma = 1 1      * gamma prior for alpha

   rgene_gamma = 2 20 1    * gammaDir prior for rate for genes

  sigma2_gamma = 1 10 1   * gammaDir prior for sigma^2     (for clock=2 or 3)

      finetune = 1: .3 .2 .6 .99 .9 .1 * auto (0 or 1): times, musigma2, rates, mixing, paras, FossilErr

         print = 1   * 0: no mcmc sample; 1: everything except branch rates 2: everything

        burnin = 20000

      sampfreq = 2

       nsample = 200000

*** Note: Make your window wider (100 columns) before running the program.

         

     

 the partial output log file was as follows:

  5% 0.36 0.14 0.13 0.43 0.43  17.635 4.282 3.961 3.393 1.472 0.897 - 1.720 0.949 -3296679.3 1:28:33

 10% 0.35 0.29 0.26 0.44 0.35  17.354 4.461 3.830 3.200 1.617 1.136 - 1.557 0.949 -3296678.8 2:12:20

 15% 0.35 0.25 0.26 0.47 0.31  17.449 4.153 3.631 2.862 1.666 1.174 - 1.463 0.949 -3296678.3 2:56:06

 20% 0.35 0.28 0.27 0.44 0.33  17.603 3.962 3.448 2.777 1.692 1.075 - 1.420 0.949 -3296678.0 3:40:05

 25% 0.35 0.32 0.29 0.47 0.31  17.613 4.010 3.337 2.552 1.582 1.044 - 1.382 0.949 -3296677.9 4:24:04

 30% 0.35 0.31 0.31 0.45 0.29  17.700 3.902 3.270 2.391 1.538 1.016 - 1.354 0.949 -3296677.8 5:08:09

 35% 0.35 0.32 0.30 0.42 0.28  17.757 3.751 3.141 2.356 1.568 1.049 - 1.341 0.949 -3296677.8 6:18:03

 40% 0.35 0.35 0.33 0.41 0.26  17.718 3.693 3.054 2.289 1.570 1.063 - 1.329 0.949 -3296677.7 7:02:32

 45% 0.35 0.38 0.34 0.41 0.27  17.679 3.605 2.931 2.203 1.534 1.026 - 1.311 0.949 -3296677.7 7:52:45

 50% 0.35 0.35 0.32 0.39 0.26  17.704 3.448 2.825 2.148 1.515 1.021 - 1.301 0.949 -3296677.6 8:43:39

 55% 0.36 0.35 0.31 0.39 0.25  17.716 3.337 2.744 2.094 1.485 1.016 - 1.288 0.949 -3296677.6 9:27:45

 60% 0.36 0.34 0.30 0.39 0.24  17.700 3.275 2.711 2.033 1.452 1.010 - 1.279 0.949 -3296677.6 10:34:43

 65% 0.36 0.33 0.31 0.38 0.23  17.670 3.224 2.679 1.989 1.444 1.021 - 1.272 0.949 -3296677.5 11:18:28

 70% 0.36 0.33 0.32 0.39 0.22  17.681 3.192 2.656 1.970 1.438 1.027 - 1.275 0.949 -3296677.5 12:54:08

 75% 0.36 0.33 0.33 0.38 0.22  17.663 3.158 2.632 1.963 1.448 1.044 - 1.273 0.949 -3296677.5 14:24:49

 80% 0.36 0.34 0.33 0.38 0.22  17.647 3.173 2.643 1.987 1.465 1.054 - 1.273 0.949 -3296677.5 15:08:34

 85% 0.36 0.34 0.33 0.37 0.21  17.643 3.139 2.617 1.977 1.471 1.068 - 1.271 0.949 -3296677.5 15:52:32

 90% 0.36 0.35 0.34 0.36 0.21  17.634 3.100 2.565 1.945 1.458 1.068 - 1.267 0.949 -3296677.4 16:36:48

 95% 0.36 0.34 0.33 0.36 0.20  17.635 3.029 2.512 1.917 1.440 1.063 - 1.263 0.949 -3296677.4 17:20:34

100% 0.36 0.33 0.33 0.36 0.20  17.671 2.967 2.469 1.893 1.424 1.053 - 1.256 0.949 -3296677.4 18:04:29

the propotion can be acceptance.

 How should  I revise my parameters to make my results reasonable? 

[Ziheng]

there is a bug in some recent versions, as described below.  Perhaps you used a version with the bug.  can you check and get the recent update 4.9e.

sorry about the trouble.

ziheng

Version 4.9e, March 2017

(*) Edited the readme files to change the license to GPL.

(*) mcmctree.  A bug was introduced in version 4.9b which causes the
program to read the fossil calibration information in the tree file
incorrectly, if joint (minimum and maximum) bounds are specified using
the symbol '<' and '>'.  If you use the notation "B()", "L()", and
'U()', the information is read correctly.  This bug was introduced in
version 4.9b and exists in 4.9c and 4.9d.  Versions 4.9a and earlier
were correct.

[Jenna Mccullough]

Hi Ziheng, 

I appear to be having this same problem. I've tried it with the ">" version and also the L() notation. But when I try to use the L() notation, I get extremely odd partial output (the calibration is for .54 hundred million years, I'm getting estimates of >2.3). I am using 4.9e. 

[Ziheng]

L means lower (minimum-age) bound.  do you have some maximum-age bounds?  you need both to push the ages to reasonable ranges.  if you have minimum bounds, the ages can be arbitrarily old.

ziheng

[Jenna Mccullough]

Ziheng, 

I have 3 calibrations: an upper bound on the first node, a node with both upper and lower bounds, and a node with a lower bound. I attached the input tree. 

I also realized that the odd numbers are because I needed to rerun the analysis with usedata=3. I have since done that, and this is the preliminary output for MCMCtree for the usedata=2 analysis. This is being run on Paml 4.8. 

I'm worried about the 5th acceptance proportion value. 

 5% 0.40 0.41 0.40 0.40 0.00  0.644 0.631 0.604 0.584 0.520 - 0.181 -189.7 12:56:03

 10% 0.40 0.41 0.40 0.40 0.00  0.645 0.632 0.603 0.583 0.520 - 0.180 -190.5 18:36:25

 15% 0.40 0.41 0.40 0.40 0.00  0.646 0.632 0.603 0.582 0.515 - 0.174 -154.4 23:57:40

 20% 0.40 0.41 0.40 0.40 0.00  0.646 0.632 0.603 0.582 0.517 - 0.177 -210.8 29:10:48

 25% 0.40 0.41 0.40 0.40 0.00  0.646 0.632 0.603 0.582 0.518 - 0.178 -180.5 34:06:47

 30% 0.40 0.41 0.40 0.40 0.00  0.646 0.633 0.603 0.582 0.517 - 0.173 -184.6 39:02:48

[Ziheng]

i think  the 5th acceptance proportion value is 0 because the move does not exist.  this is a move that changes parameters in the substitution model and because you are using approximate likelihood calculation, the move is not used.

please simply ignore that number.

ziheng