etetoolkit / ete

Python package for building, comparing, annotating, manipulating and visualising trees. It provides a comprehensive API and a collection of command line tools, including utilities to work with the NCBI taxonomy tree.
http://etetoolkit.org
GNU General Public License v3.0
782 stars 212 forks source link

Availability of runmode = -3 in codeml? [bayesian estimate of w] #198

Open BitaoQiu opened 8 years ago

BitaoQiu commented 8 years ago

Is it possible to integrate run mode = -3 in ete? This is corresponding to bayesian estimation of dn/ds in yang's paper. http://abacus.gene.ucl.ac.uk/software/pamlDOC.pdf (page 34) and according to Yang, it's more robust and reduce extreme dn/ds. "Bayesian Estimation of Nonsynonymous/Synonymous Rate Ratios for Pairwise Sequence Comparisons"

fransua commented 8 years ago

Hi @StanQiu , a priori any option can be passed to the ete3 evol tool. In your case you could try:

ete3 evol --codeml_param runmode,-3 ...

note that this option will be applied to all the models to be run.

Note that I am assuming you are using the command line tool, in case I am wrong you can set the same parameter through the run_model function (http://etetoolkit.org/docs/latest/reference/reference_evoltree.html#ete3.EvolNode.run_model)

let us know if it helps

BitaoQiu commented 8 years ago

Hi @fransua,

Thanks for replying. I tried but failed...and it says something wrong with the codeml_parser. I checked the code of ete3 and it seems the codeml parser is designed for run mode = -2 (ML)

Adding some comments...

Sorry for the trouble. I read the paper of Yang again and realize that run mode = -3 will only be useful when doing pairwise dn/ds calculation, where the result is already available in the rst file. So even ete returns error, it's already in the final step and will not influence the analysis...

Running CodeML/Slr (0 CPUs)
  - processing model M0
lnL =-1354.183487
test/M0/out
Traceback (most recent call last):
  File "/seqdata/krogh/cse/bitao/bin/anaconda_ete/bin/ete3", line 9, in <module>
    load_entry_point('ete3==3.0.0b33', 'console_scripts', 'ete3')()
  File "/seqdata/krogh/cse/bitao/bin/anaconda_ete/lib/python2.7/site-packages/ete3/tools/ete.py", line 92, in main
    _main(sys.argv)
  File "/seqdata/krogh/cse/bitao/bin/anaconda_ete/lib/python2.7/site-packages/ete3/tools/ete.py", line 238, in _main
    args.func(args)
  File "/seqdata/krogh/cse/bitao/bin/anaconda_ete/lib/python2.7/site-packages/ete3/tools/ete_evol.py", line 874, in run
    load_model(model, tree, models[model], **params)
  File "/seqdata/krogh/cse/bitao/bin/anaconda_ete/lib/python2.7/site-packages/ete3/tools/ete_evol.py", line 667, in load_model
    tree.link_to_evol_model(path, model_obj)
  File "/seqdata/krogh/cse/bitao/bin/anaconda_ete/lib/python2.7/site-packages/ete3/evol/evoltree.py", line 421, in link_to_evol_model
    model._load(path)
  File "/seqdata/krogh/cse/bitao/bin/anaconda_ete/lib/python2.7/site-packages/ete3/evol/model.py", line 163, in _load
    parse_paml(path, self)
  File "/seqdata/krogh/cse/bitao/bin/anaconda_ete/lib/python2.7/site-packages/ete3/evol/parser/codemlparser.py", line 249, in parse_paml
    model.stats ['np' ] = int (float  (line.split()[0]))
ValueError: could not convert string to float: lnL
fransua commented 8 years ago

Ok, I see. Yes, the module is not prepared for this feature. We will have to think about it (what to keep and display), even if it's recommended for paiwise dn/ds, the feature and the results of the runmode=-3 are interesting. thanks