Bug in rootAt

[Bastien Boussau]

Hi all,

I found a bug in rootAt. Basically when using rootAt, a branch length (leading to the node used as outgroup) disappears. This function used to work for me: it was used in an old piece of code, that became non-functional.

I think the problem lies in line 501 in TreeTemplate.h:
newRoot->deleteDistanceToFather();

If I remove it, my code functions again. I do not understand why this line was necessary, but I have not given it a huge amount of thought.

I tried to commit this change but failed:
fatal: Unable to create '/git/bpp-phyl.git/gc.pid.lock': Permission non accordée
To ssh://bou...@biopp.univ-montp2.fr/git/bpp-phyl.git/
9a398b4..b64f87a master -> master

I may be using the wrong server or something.

If someone wants to commit this change for me, that’s great!

Thanks,

Best,

Bastien

[Julien Yann Dutheil]

Hi Bastien,

I actually added this line recently to solve another bug... I did not see this side effect, sorry! The deletedistancetofather has to be here. Will fix rootat accordingly asap.

Julien.

--
You received this message because you are subscribed to the Google Groups "Bio++ Development Forum" group.
To unsubscribe from this group and stop receiving emails from it, send an email to biopp-devel-fo...@googlegroups.com.
To post to this group, send email to biopp-de...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/biopp-devel-forum/8653EDE5-1D57-4A55-AE5A-7CA3219FE1C7%40gmail.com.
For more options, visit https://groups.google.com/d/optout.

[Bastien Boussau]

Great, thanks very much Julien. Good thing I was not able to push my change, it would have been annoying for you.
Bastien

To view this discussion on the web visit https://groups.google.com/d/msgid/biopp-devel-forum/CAFr7pNUxp%2BJ9pKkLAM172d-T2AG3qHcDRR0%2B7Q%2B%3DAnZW3_x5Mw%40mail.gmail.com.

[Julien Yann Dutheil]

Hi Bastien,

First of all, your commit worked i think. Second, unlike what i thought, this is not related to the other bug i fixed. Will need to dig a bit more in the git history to investigate this... I keep you updated.

J.

To view this discussion on the web visit https://groups.google.com/d/msgid/biopp-devel-forum/CABHONyxunxFoA7qk1h97D5HZ9D%2B6PiCCBW0S39irW4f%2B3jO6wg%40mail.gmail.com.

[Bastien Boussau]

Ok, thanks, sorry it's turning out to be more complicated...
Bastien

To view this discussion on the web visit https://groups.google.com/d/msgid/biopp-devel-forum/CAFr7pNUsxpUThOUWeRqV9w9Q_RCX%2BogshWmog7xHMzh3t8RVdw%40mail.gmail.com.

[Julien Yann Dutheil]

Hi Bastien,

That is strange, indeed this line should not be necessary, as root nodes are not supposed to have branch lengths. Are you sure the issue come from this function :s ? I can see it uses the getPathBetweenAnyTwoNodes method that we recently updated, so that could be a potential source of pb. I will write a small unit test for the rootAt function and see if that helps.

J.

To view this discussion on the web visit https://groups.google.com/d/msgid/biopp-devel-forum/CABHONyzLQ2saLxKnYfKVy1FLuDNZ_xbGwKN425dL5_sObLctAg%40mail.gmail.com.

For more options, visit https://groups.google.com/d/optout.

--

Julien Y. Dutheil, Ph-D

0 (+49) 4522 763 298

§ Max Planck Institute for Evolutionary Biology

Molecular Systems Evolution

Department of Evolutionary Genetics

Plön -- GERMANY

§ Intitute of Evolutionary Sciences - Montpellier
University of Montpellier 2 -- FRANCE

[Bastien Boussau]

Hi Julien,

Here is the function rootAt with my “prints” in it: 

void rootAt(N* newRoot)

  {

    if (root_ == newRoot) return;

    if (isRooted()) unroot();

    std::vector<Node*> path = TreeTemplateTools::getPathBetweenAnyTwoNodes(*root_, *newRoot);

    std::cout << "rootAt: " << TreeTemplateTools::treeToParenthesis(*this, true) <<std::endl;

    std::cout << "newRoot->getName(): "<< newRoot->getName() <<std::endl;

    for (size_t i = 0; i < path.size() - 1; i++)

    {

      // pathMatrix[i] -> _father = pathMatrix[i + 1];

      // pathMatrix[i] -> setDistanceToFather(pathMatrix[i + 1] -> getDistanceToFather());

      // typename vector<Node *>::iterator vec_iter;

      // vec_iter = remove(pathMatrix[i] -> _sons.begin(), pathMatrix[i] -> _sons.end(), pathMatrix[i + 1]);

      // pathMatrix[i] -> _sons.erase(vec_iter, pathMatrix[i] -> _sons.end()); // pg 1170, primer.

      // pathMatrix[i+1] -> _sons.push_back(pathMatrix[i + 1] -> getFather());

      // pathMatrix[i+1] -> _father = 0;

      path[i]->removeSon(path[i + 1]);

//std::cout << " dist to father "<< path[i]->getDistanceToFather() << std::endl;

      if (path[i + 1]->hasDistanceToFather())  { 

std::cout << "HEHE " <<std::endl;

path[i]->setDistanceToFather(path[i + 1]->getDistanceToFather());

}

      else path[i]->deleteDistanceToFather();

      path[i + 1]->addSon(path[i]);

      std::vector<std::string> names = path[i + 1]->getBranchPropertyNames();

      for (size_t j = 0; j < names.size(); j++)

      {

        path[i]->setBranchProperty(names[j], *path[i + 1]->getBranchProperty(names[j]));

      }

      path[i + 1]->deleteBranchProperties();

    }

std::cout << " dist to father2 "<< newRoot->getDistanceToFather() << std::endl;

    newRoot->deleteDistanceToFather();

std::cout << " dist to father3 "<< newRoot->getDistanceToFather() << std::endl;

    newRoot->deleteBranchProperties();

    root_ = newRoot;

  }

And here is the output I got: 

rootAt: (a1_0:0.01111,(a2_1:0.01,a3_2:0.01)3:0.01,((b1_5:0.01,(b2_6:0.01,b3_7:0.01)8:0.01)9:0.01,(c1_10:0.01,(c2_11:0.01,c3_12:0.01)13:0.01)14:0.01)15:0.02);

newRoot->getName(): a1

HEHE 

 dist to father2 0.01111

_____________________________________________________

ERROR!!!

NodeException: Node::getDistanceToFather: Node has no distance.(id:0)

The starting tree was this one: 

((a1:0.01111,(a2:0.01,a3:0.01):0.01):0.01,((b1:0.01,(b2:0.01,b3:0.01):0.01):0.01,(c1:0.01,(c2:0.01,c3:0.01):0.01):0.01):0.01);

That’s all the evidence I have to base my claim that the problem was coming from the "newRoot->deleteDistanceToFather();” before  std::cout << " dist to father3 "<< newRoot->getDistanceToFather() << std::endl;

But I must confess I’m not sure what’s happening in this function.

Bastien

To view this discussion on the web visit https://groups.google.com/d/msgid/biopp-devel-forum/CAFr7pNWmTLWxS72jT2RUtPSOpbJpexYh7MmWScdAT9moPnVijQ%40mail.gmail.com.

[Julien Yann Dutheil]

Hi Bastien,

I have updated the corresponding code a bit. I did not find any major issue there though, so I cannot warranty this would solve your probelm as well :s

I have created a unit test for the rootAt function. It creates a random tree and compute the total branch length. Then it pick a node randomly (a 100 times), change the root, and recompute the total length. If it is not equal to the initial sum, the test fails. So far I could not see any misbehavior...

Please keep me updated :)

Cheers,

Julien.