1/*************************************************************************************/ 2 3function PostOrderAVL2StringDL (theAVL, doLengths) 4{ 5 return PostOrderAVL2StringAnnotate (theAVL, doLengths, ""); 6} 7 8/*************************************************************************************/ 9 10function PostOrderAVL2StringAnnotate (theAVL, doLengths,label) 11{ 12 return PostOrderAVL2StringAnnotateAux (theAVL, doLengths, label, "[]"); 13} 14 15/*************************************************************************************/ 16 17function PostOrderAVL2StringAnnotateAux (theAVL, doLengths, label, chars) 18{ 19 _ost = ""; 20 _ost * 256; 21 22 lastLevel = 0; 23 treeSize = Abs(theAVL); 24 treeInfo = theAVL[0]; 25 rootIndex = treeInfo["Root"]; 26 27 for (nodeIndex = 1; nodeIndex < treeSize; nodeIndex = nodeIndex + 1) 28 { 29 nodeInfo = theAVL[nodeIndex]; 30 myDepth = nodeInfo["Depth"]; 31 if (lastDepth < myDepth) 32 { 33 if (lastDepth) 34 { 35 _ost * ","; 36 } 37 for (pidx = lastDepth; pidx < myDepth; pidx = pidx + 1) 38 { 39 _ost * "("; 40 } 41 } 42 else 43 { 44 if (lastDepth > myDepth) 45 { 46 for (pidx = myDepth; pidx < lastDepth; pidx = pidx + 1) 47 { 48 _ost * ")"; 49 } 50 } 51 else 52 { 53 _ost * ","; 54 } 55 } 56 57 _ost * nodeInfo["Name"]; 58 59 if (Abs (label)) 60 { 61 if (Abs(nodeInfo[label])) 62 { 63 _ost * (chars[0] + nodeInfo[label] + chars[1]); 64 } 65 } 66 67 if (doLengths) 68 { 69 if (nodeIndex < treeSize - 1) 70 { 71 _ost * ":"; 72 _ost * (""+nodeInfo ["Length"]); 73 } 74 } 75 lastDepth = myDepth; 76 } 77 78 _ost * 0; 79 return _ost; 80} 81 82 83/*************************************************************************************/ 84 85function PostOrderAVL2String (theAVL) 86{ 87 return PostOrderAVL2StringDL(theAVL, 0); 88} 89 90/*************************************************************************************/ 91 92function PostOrderAVL2StringDistances (theAVL, distAVL) 93{ 94 _ost = ""; 95 _ost * 256; 96 97 lastLevel = 0; 98 treeSize = Abs(theAVL); 99 treeInfo = theAVL[0]; 100 rootIndex = treeInfo["Root"]; 101 102 for (nodeIndex = 1; nodeIndex < treeSize; nodeIndex+=1) 103 { 104 nodeInfo = theAVL[nodeIndex]; 105 myDepth = nodeInfo["Depth"]; 106 myName = nodeInfo["Name"]; 107 108 if (lastDepth < myDepth) 109 { 110 if (lastDepth) 111 { 112 _ost * ","; 113 } 114 for (pidx = lastDepth; pidx < myDepth; pidx = pidx + 1) 115 { 116 _ost * "("; 117 } 118 } 119 else 120 { 121 if (lastDepth > myDepth) 122 { 123 for (pidx = myDepth; pidx < lastDepth; pidx = pidx + 1) 124 { 125 _ost * ")"; 126 } 127 } 128 else 129 { 130 _ost * ","; 131 } 132 } 133 if (nodeIndex != rootIndex) { 134 _ost * myName; 135 } 136 if (nodeIndex < treeSize - 1) 137 { 138 _ost * ":"; 139 _ost * (""+distAVL [myName]); 140 } 141 lastDepth = myDepth; 142 } 143 144 _ost * 0; 145 return _ost; 146} 147 148/*************************************************************************************/ 149 150function KillInternalZeroBranchLengths (treeAVL) 151{ 152 treeSize = Abs(treeAVL); 153 if (treeSize == 3) 154 { 155 return "(" + (treeAVL[0])["Name"] + "," + (treeAVL[1])["Name"] + ")"; 156 } 157 newTreeAVL = {}; 158 oldIndexMap= {treeSize,1}; 159 index2 = {treeSize,1}; 160 newDAVL = {}; 161 newTreeAVL [0] = treeAVL[0]; 162 allDeleted = 0; 163 for (nodeIndex = 1; nodeIndex < treeSize; nodeIndex = nodeIndex + 1) 164 { 165 newDAVL [(treeAVL[nodeIndex])["Name"]] = (treeAVL[nodeIndex])["Length"]; 166 167 if (Abs((treeAVL[nodeIndex])["Children"]) && Abs((treeAVL[nodeIndex])["Length"]) <= 1e-10 && (treeAVL[nodeIndex])["Parent"]) 168 /* zero internal branch */ 169 { 170 oldIndexMap[nodeIndex] = -(treeAVL[nodeIndex])["Parent"]; 171 allDeleted = allDeleted + 1; 172 } 173 else 174 { 175 newTreeAVL [nodeIndex-allDeleted] = treeAVL[nodeIndex]; 176 oldIndexMap[nodeIndex] = nodeIndex-allDeleted; 177 index2 [nodeIndex-allDeleted] = nodeIndex; 178 } 179 } 180 181 if (allDeleted) 182 { 183 markedIndices = {}; 184 for (nodeIndex = treeSize-1; nodeIndex>0; nodeIndex = nodeIndex - 1) 185 { 186 if (oldIndexMap [nodeIndex]<0) 187 { 188 markedIndices[nodeIndex] = 1; 189 oldIndexMap[nodeIndex] = oldIndexMap[-oldIndexMap[nodeIndex]]; 190 } 191 } 192 treeSize = Abs (newTreeAVL); 193 for (nodeIndex = 1; nodeIndex<treeSize; nodeIndex = nodeIndex + 1) 194 { 195 meParent = (newTreeAVL[nodeIndex])["Parent"]; 196 _cc = Abs((newTreeAVL[nodeIndex])["Children"]); 197 if (_cc > 0) 198 { 199 newChildrenMap = {}; 200 for (_cci = 0; _cci < _cc; _cci = _cci+1) 201 /* map children to new indices */ 202 { 203 _cn = ((newTreeAVL[nodeIndex])["Children"])[_cci]; 204 if (markedIndices[_cn] == 0) 205 { 206 newChildrenMap[Abs(newChildrenMap)] = oldIndexMap[_cn]; 207 } 208 } 209 ((newTreeAVL[nodeIndex])["Children"]) = newChildrenMap; 210 } 211 212 if (meParent > 0) 213 { 214 meParentOI = meParent; 215 meParent = oldIndexMap[meParent]; 216 (newTreeAVL[nodeIndex])["Parent"] = meParent; 217 if (markedIndices[meParentOI]) 218 { 219 /* 220 fprintf (stdout, "Insert ", (newTreeAVL[nodeIndex])["Name"], " as a child of ", (newTreeAVL[meParent])["Name"], " index ", index2[nodeIndex], "(", oldIndexMap[index2[nodeIndex]], ",", nodeIndex,")\n"); 221 */ 222 ((newTreeAVL[meParent])["Children"])[Abs((newTreeAVL[meParent])["Children"])] = index2[nodeIndex]; 223 } 224 } 225 226 } 227 228 229 for (nodeIndex = treeSize-1; nodeIndex>0; nodeIndex = nodeIndex - 1) 230 { 231 _cc = Abs((newTreeAVL[nodeIndex])["Children"]); 232 if (_cc > 0) 233 { 234 _cd = (newTreeAVL[nodeIndex])["Depth"] + 1; 235 for (_cci = 0; _cci < _cc; _cci = _cci+1) 236 { 237 _cn = ((newTreeAVL[nodeIndex])["Children"])[_cci]; 238 /*fprintf (stdout, (newTreeAVL[_cn])["Name"], ":", (newTreeAVL[_cn])["Depth"], "=>", _cd, "\n");*/ 239 (newTreeAVL[_cn])["Depth"] = _cd; 240 } 241 } 242 } 243 244 (newTreeAVL[0])["Root"] = treeSize-1; 245 } 246 247 return PostOrderAVL2StringDistances (newTreeAVL, newDAVL); 248} 249 250 251/*************************************************************************************/ 252 253function TreeAVL2String (treeAVL) 254{ 255 rootNode = treeAVL[0]; 256 rootNode = rootNode["Root"]; 257 return subtreeAVLStr (rootNode,0,0); 258} 259 260 261/*************************************************************************************/ 262 263function subtreeAVLStr (nodeIndex,k,treeString) 264{ 265 nodeInfo = treeAVL[nodeIndex]; 266 k = Abs(nodeInfo["Children"])-1; 267 if (k>=0) 268 { 269 while (k>=0) 270 { 271 nodeInfo = treeAVL[nodeIndex]; 272 cNodes = nodeInfo["Children"]; 273 cNodes = cNodes[k]; 274 if (k < Abs(nodeInfo["Children"])-1) 275 { 276 ExecuteCommands("treeString=subtreeAVLStr (cNodes,k,treeString)+\",\"+treeString;"); 277 } 278 else 279 { 280 ExecuteCommands("treeString=subtreeAVLStr (cNodes,k,treeString)+\")\";"); 281 } 282 k=k-1; 283 } 284 return "("+treeString; 285 } 286 else 287 { 288 callLevel = callLevel - 1; 289 return nodeInfo["Name"]; 290 } 291} 292 293/*************************************************************************************/ 294 295function InsertANode (theAVL&,insertAt,newNodeName) 296{ 297 nodeInfo = theAVL[insertAt]; 298 if (Abs(nodeInfo)) 299 { 300 nparent = nodeInfo["Parent"]; 301 if (nparent > 0) 302 { 303 lastIndex = Abs(theAVL); 304 myDepth = nodeInfo["Depth"]; 305 newParentNode = {}; 306 newParentNode ["Name"] = "Node"+lastIndex; 307 newParentNode ["Parent"] = nparent; 308 newParentNode ["Depth"] = myDepth; 309 310 newChildNode = {}; 311 newChildNode ["Name"] = newNodeName; 312 newChildNode ["Parent"] = lastIndex; 313 newChildNode ["Depth"] = myDepth + 1; 314 315 pChildren = {}; 316 pChildren [0] = insertAt; 317 pChildren [1] = lastIndex+1; 318 newParentNode ["Children"] = pChildren; 319 320 theAVL[lastIndex] = newParentNode; 321 theAVL[lastIndex+1] = newChildNode; 322 323 /* update the parent*/ 324 325 nodeInfo ["Parent"] = lastIndex; 326 theAVL[insertAt] = nodeInfo; 327 328 /* update the list of children at the parent node*/ 329 330 parentInfo = theAVL[nparent]; 331 parentChildren = parentInfo["Children"]; 332 333 for (nic = Abs(parentChildren)-1; nic >= 0; nic = nic-1) 334 { 335 if (parentChildren[nic] == insertAt) 336 { 337 break; 338 } 339 } 340 341 parentChildren[nic] = lastIndex; 342 parentInfo["Children"] = parentChildren; 343 theAVL[nparent] = parentInfo; 344 345 /* now update the depths at new NodeName and all of its children */ 346 347 nodeCache = {}; 348 nodeCache[0] = insertAt; 349 cacheIndex = 0; 350 351 while (cacheIndex <= Abs(nodeCache)) 352 { 353 nparent = nodeCache[cacheIndex]; 354 nodeInfo = theAVL[nparent]; 355 nodeInfo["Depth"] = nodeInfo["Depth"] + 1; 356 theAVL[nparent] = nodeInfo; 357 nodeChildren = nodeInfo["Children"]; 358 for (nic = Abs(nodeChildren)-1; nic >=0; nic = nic-1) 359 { 360 nodeCache [Abs(nodeCache)] = nodeChildren[nic]; 361 } 362 cacheIndex = cacheIndex + 1; 363 } 364 365 nodeCache = 0; 366 } 367 } 368 return 0; 369} 370 371/*************************************************************************************/ 372 373function ModifyDepth (nIndex, modAmount) 374{ 375 nodeInfo = theAVL[nIndex]; 376 nodeInfo ["Depth"] = nodeInfo ["Depth"] + modAmount; 377 theAVL[nIndex] = nodeInfo; 378 379} 380 381/*************************************************************************************/ 382 383function echoAVL (anAVL) 384{ 385 for (k=1; k<Abs(anAVL); k=k+1) 386 { 387 nodeInfo = anAVL[k]; 388 myChildren = nodeInfo["Children"]; 389 if (Abs(myChildren)) 390 { 391 fprintf (stdout, "Index ", k, ":", nodeInfo["Name"], " : parent = ", nodeInfo["Parent"], " children:"); 392 for (k2 = 0; k2 < Abs(myChildren); k2=k2+1) 393 { 394 fprintf (stdout,"\t", myChildren[k2]); 395 } 396 fprintf (stdout, " depth: ", nodeInfo["Depth"], "\n"); 397 398 } 399 else 400 { 401 fprintf (stdout, "Index ", k, ":", nodeInfo["Name"], " : parent = ", nodeInfo["Parent"], " children: none, depth: ", nodeInfo["Depth"], "\n"); 402 } 403 } 404 return 0; 405} 406 407/*************************************************************************************/ 408 409function selectATreeBranch (treeID, title) 410{ 411 ExecuteCommands ("internalNodes = BranchCount(`treeID`);"); 412 ExecuteCommands ("leafNodes = TipCount(`treeID`);"); 413 414 choiceMatrix = {internalNodes+leafNodes,2}; 415 416 for (bc=0; bc<internalNodes; bc=bc+1) 417 { 418 ExecuteCommands ("choiceMatrix[bc][0] = BranchName(`treeID`,bc);choiceMatrix[bc][1] = \"Internal Branch Rooting \" + `treeID`[bc];"); 419 } 420 421 for (bc=0; bc<leafNodes; bc=bc+1) 422 { 423 ExecuteCommands ("choiceMatrix[bc+internalNodes][0] = TipName(`treeID`,bc);"); 424 choiceMatrix[bc+internalNodes][1] = "Terminal branch endin in " + choiceMatrix[bc+internalNodes][0]; 425 } 426 427 mxTreeSpec = {5,1}; 428 mxTreeSpec [0] = treeID; 429 mxTreeSpec [1] = "8240"; 430 mxTreeSpec [2] = "10,40,-10,-175,1"; 431 mxTreeSpec [3] = ""; 432 mxTreeSpec [4] = ""; 433 OpenWindow (TREEWINDOW, mxTreeSpec); 434 435 ChoiceList (bOption,title,1,NO_SKIP,choiceMatrix); 436 437 if (bOption < 0) 438 { 439 return ""; 440 } 441 return choiceMatrix[bOption][0]; 442} 443 444/*************************************************************************************/ 445 446function computeTreeSplits (treeID,mirror,leafMap) 447{ 448 ExecuteCommands ("_treeAVL = `treeID`^0;"); 449 ExecuteCommands ("_leafNodes = TipCount(`treeID`);"); 450 _tipMap = {}; 451 _splitMap = {}; 452 _splitTemplate = {_leafNodes,1}; 453 _hasLeafMap = Abs (leafMap); 454 455 for (_k = 1; _k < Abs(_treeAVL); _k = _k+1) 456 { 457 if (Abs((_treeAVL[_k])["Children"]) == 0) 458 { 459 if (_hasLeafMap) 460 { 461 _tipMap [_k] = leafMap [(_treeAVL[_k])["Name"]]; 462 } 463 else 464 { 465 _tipMap [_k] = Abs(_tipMap); 466 } 467 } 468 (_treeAVL[_k])["Split"] = _splitTemplate; 469 } 470 471 472 for (_k = 1; _k < Abs(_treeAVL); _k = _k+1) 473 { 474 _pc = (_treeAVL[_k])["Parent"]; 475 if (_pc) 476 { 477 _cc = Abs((_treeAVL[_k])["Children"]); 478 if (_cc) 479 { 480 _mySplit = ((_treeAVL[_k])["Split"])["_MATRIX_ELEMENT_VALUE_>0"]; 481 (_treeAVL[_pc])["Split"] = (_treeAVL[_pc])["Split"] + _mySplit; 482 } 483 else 484 { 485 ((_treeAVL[_k])["Split"])[_tipMap[_k]] = 1; 486 _mySplit = (_treeAVL[_k])["Split"]; 487 ((_treeAVL[_pc])["Split"])[_tipMap[_k]] = 1; 488 } 489 _splitMap[_mySplit] = (_treeAVL[_k])["Name"]; 490 } 491 } 492 493 494 _uniqueSplits = Rows(_splitMap); 495 _splitStrings = {}; 496 for (_k = 0; _k < Abs(_splitMap); _k = _k+1) 497 { 498 _stringKey = ""; _stringKey * 128; 499 _stringKey2 = ""; _stringKey2 * 128; 500 ExecuteCommands ("_thisKey = " + _uniqueSplits[_k]); 501 for (_k2 = 0; _k2 < _leafNodes; _k2 = _k2 + 1) 502 { 503 if (_thisKey[_k2]) 504 { 505 _stringKey * "*"; 506 _stringKey2 * "-"; 507 } 508 else 509 { 510 _stringKey * "-"; 511 _stringKey2 * "*"; 512 } 513 } 514 _stringKey * 0;_stringKey2 * 0; 515 _splitStrings [_stringKey] = _splitMap[_uniqueSplits[_k]]; 516 if (mirror) 517 { 518 _splitStrings [_stringKey2] = _splitStrings [_stringKey]; 519 } 520 } 521 DeleteObject (_splitMap); 522 return _splitStrings; 523} 524 525/*************************************************************************************/ 526 527/* 528 return the most recent common ancestor for a group 529 of nodes in an AVL; returns the index of the MRCA 530*/ 531 532function _findMRCA (treeAVL, nodeIDList) 533{ 534 _nodeCount = Rows(nodeIDList)*Columns(nodeIDList); 535 536 if (_nodeCount) 537 { 538 nodeIDList = nodeIDList % 0; 539 _highestNode = (treeAVL[nodeIDList[0]])["Depth"]; 540 _currentDepth = _highestNode; 541 _lastNode = nodeIDList[_nodeCount-1]; 542 for (_nodeIterator = nodeIDList[0]+1; _nodeIterator <= _lastNode; 543 _nodeIterator = _nodeIterator + 1) 544 { 545 _currentDepth = (treeAVL[_nodeIterator])["Depth"]; 546 if (_currentDepth < _highestNode) 547 { 548 _highestNode = _currentDepth; 549 } 550 } 551 552 for (; _nodeIterator <= Abs(treeAVL); _nodeIterator = _nodeIterator+1) 553 { 554 if ((treeAVL[_nodeIterator])["Depth"] < _highestNode) 555 { 556 return _nodeIterator; 557 } 558 } 559 560 } 561 562 return Abs(treeAVL)-1; 563} 564 565/*************************************************************************************/ 566/* 567 use parsimony to reconstruct ancestral states based on leaf AVL labels 568 return the number of substitutions 569*/ 570 571function _parsimonyAncestralMapping (treeAVL&, label) 572{ 573 uniqueLabels = {}; 574 idToLabel = {}; 575 for (_nodeIterator = 1; _nodeIterator < Abs (treeAVL); _nodeIterator = _nodeIterator + 1) 576 { 577 isLeaf = Abs((treeAVL[_nodeIterator])["Children"]) == 0; 578 if (isLeaf) 579 { 580 isLeaf = (treeAVL[_nodeIterator])[label]; 581 if (uniqueLabels[isLeaf] == 0) 582 { 583 uniqueLabels [isLeaf] = 1+Abs(uniqueLabels); 584 idToLabel[Abs(idToLabel)] = isLeaf; 585 } 586 } 587 } 588 589 assignmentMatrices = {}; 590 kindCount = Abs (uniqueLabels); 591 592 for (_nodeIterator = 1; _nodeIterator < Abs (treeAVL); _nodeIterator = _nodeIterator + 1) 593 { 594 aMx = {2,kindCount}; 595 s2 = uniqueLabels[(treeAVL[_nodeIterator])[label]]-1; 596 for (k=0; k<kindCount; k=k+1) 597 { 598 aMx[0][k] = s2; 599 aMx[1][k] = 1-(k==s2); 600 } 601 assignmentMatrices [(treeAVL[_nodeIterator])["Name"]] = aMx; 602 } 603 604 for (_nodeIterator = 1; _nodeIterator < Abs (treeAVL); _nodeIterator = _nodeIterator + 1) 605 { 606 nodeInfo = treeAVL[_nodeIterator]; 607 nodeChildren = nodeInfo ["Children"]; 608 cCount = Abs(nodeChildren); 609 610 if (cCount) 611 { 612 localMatrices = {}; 613 614 nodeName = nodeInfo["Name"]; 615 616 for (s1 = 0; s1<cCount; s1=s1+1) 617 { 618 localMatrices[s1] = assignmentMatrices[(treeAVL [nodeChildren[s1]]) ["Name"]]; 619 } 620 621 twoWay = {kindCount,1}; 622 623 for (s2 = 0; s2 < kindCount; s2 = s2+1) 624 { 625 lc = 0; 626 for (s3 = 0; s3<cCount; s3=s3+1) 627 { 628 lc = lc + (localMatrices[s3])[1][s2]; 629 } 630 twoWay[s2] = lc; 631 } 632 633 if (nodeInfo["Parent"]) 634 { 635 aMx = {2,kindCount}; 636 637 for (s2 = 0; s2 < kindCount; s2 = s2+1) 638 { 639 minV = 1e100; 640 minI = 0; 641 642 for (s3 = 0; s3 < kindCount; s3 = s3+1) 643 { 644 aCost = (s3!=s2) + twoWay[s3]; 645 if (minV > aCost) 646 { 647 minV = aCost; 648 minI = s3; 649 } 650 } 651 652 aMx[0][s2] = minI; 653 aMx[1][s2] = minV; 654 } 655 656 assignmentMatrices [nodeName] = aMx; 657 } 658 else 659 { 660 totalCost = 1e100; 661 rootState = 0; 662 663 for (s2 = 0; s2 < kindCount; s2 = s2+1) 664 { 665 if (twoWay[s2] < totalCost) 666 { 667 totalCost = twoWay[s2]; 668 rootState = s2; 669 } 670 } 671 (treeAVL[_nodeIterator])[label] = idToLabel[rootState]; 672 } 673 } 674 } 675 676 for (_nodeIterator = Abs (treeAVL)-1; _nodeIterator >=1 ; _nodeIterator = _nodeIterator - 1) 677 { 678 nodeInfo = treeAVL[_nodeIterator]; 679 nodeChildren = nodeInfo ["Children"]; 680 681 if (Abs(nodeChildren)) 682 { 683 nodeName = nodeInfo["Name"]; 684 nodeParent = nodeInfo["Parent"]; 685 if (nodeParent) 686 { 687 aMx = assignmentMatrices [nodeName]; 688 (treeAVL[_nodeIterator])[label] = idToLabel[aMx[0][uniqueLabels[(treeAVL[nodeParent])[label]]-1]]; 689 } 690 } 691 } 692 693 assignmentMatrices = 0; 694 return totalCost; 695} 696 697 698 699/*************************************************************************************/ 700/* 701 Center of Tree (COT) calculation 702*/ 703 704function ComputeCOT (_treeID, _filePath) 705/* IN : the identifier of an existing tree variable (String) 706 _filePath : if not an empty string, write a PostScript image with this tree 707 to the given path 708 709 OUT : an associative array with four entries 710 "Branch" - the branch where the COT resides 711 branch ENDS at the node whose value is returned 712 713 "Split" - how far up the branch the COT is, measured in the same units as branch lengths 714 from the END of the branch (i.e node) 715 716 "COTTree" - tree rerooted at the COT 717 718 "Distances" - mean squared distance from the COT to all leaves 719*/ 720{ 721 _power = 2; 722 /* '2' is the power of the distance function to 723 minimize - i.e. least squares in this case */ 724 ExecuteCommands ( 725 "_cot = Min ("+_treeID+",_power);"); 726 727 _returnList = {}; 728 _returnList ["Branch"] = _cot["COT_NODE"]; 729 _returnList ["Split"] = _cot["COT_SPLIT"]; 730 _meanD = 0; 731 _keys = Rows(_cot["COT_TO_NODE"]); 732 /* _cot["COT_TO_NODE"] is an associative array 733 mapping tree nodes to the total distance from COT (linear distance) 734 */ 735 736 for (_k = 0; _k < Columns(_keys); _k=_k+1) 737 { 738 _meanD = _meanD + (_cot["COT_TO_NODE"])[_keys[_k]]^_power; 739 } 740 741 _returnList ["Distances"] = _meanD/Columns(_keys); 742 /* reroot the tree at the COT branch */ 743 ExecuteCommands ("_rerootedTree = RerootTree ("+_treeID+", \""+_returnList ["Branch"]+"\");"); 744 /* hackish lines below split the root branch into appropriate bits */ 745 746 ACCEPT_ROOTED_TREES = 1; 747 UseModel (USE_NO_MODEL); 748 Tree _temp = _rerootedTree; 749 _tempA = _temp^0; /* this converts the tree into a post-order list of nodes as an associative array; 750 print it to see the structure */ 751 _rootID = (_tempA[0])["Root"]; 752 _dist = {}; 753 for (_k = 1; _k < Abs (_tempA); _k = _k + 1) 754 { 755 if ((_tempA[_k])["Parent"] == _rootID) 756 { 757 if (_k == _rootID - 1) 758 { 759 _dist [(_tempA[_k])["Name"]] = _cot["COT_SPLIT"]; 760 } 761 else 762 { 763 _dist [(_tempA[_k])["Name"]] = _cot["COT_BRANCH_LENGTH"]-_cot["COT_SPLIT"]; 764 } 765 } 766 else 767 { 768 _dist [(_tempA[_k])["Name"]] = (_tempA[_k])["Length"]; 769 } 770 } 771 772 _returnList["COTTree"] = PostOrderAVL2StringDistances (_tempA,_dist); 773 774 /* make a postscript file if needed */ 775 if (Abs(_filePath)) 776 { 777 TREE_OUTPUT_OPTIONS = {}; 778 nodeSpec = {}; 779 nodeSpec ["TREE_OUTPUT_BRANCH_SPLIT"] = _cot["COT_SPLIT"]/_cot["COT_BRANCH_LENGTH"]; 780 TREE_OUTPUT_OPTIONS [_cot["COT_NODE"]] = nodeSpec; 781 ExecuteCommands ("psString = PSTreeString ("+_treeID+",\"STRING_SUPPLIED_LENGTHS\",{{-1,-1}});"); 782 fprintf (_filePath, CLEAR_FILE, psString); 783 } 784 785 return _returnList; 786 787} 788 789 790/*************************************************************************************/ 791/* 792 Computes the distance from each node to the root of the tree. 793 By default, uses the "Length" key, but this can be redefined by passing an argument 794*/ 795 796function PathDistanceToRoot (_treeAVL, _distanceKey) 797/* IN : the identifier of a post-order tree AVL (dictionary) 798 _distanceKey : if not an empty string, use this key to retrieve branch lengths 799 800 OUT : an associative array keyed on branch names with paths stored as values 801 802*/ 803{ 804 _distances = {}; 805 806 if (Abs (_distanceKey) == 0) 807 { 808 _distanceKey = "Length"; 809 } 810 811 for (_k = Abs (_treeAVL)-2; _k > 0; _k = _k-1) 812 { 813 _nodeName = (_treeAVL[_k])["Name"]; 814 _distances[_nodeName] = _distances[(_treeAVL[(_treeAVL[_k])["Parent"]])["Name"]] + (_treeAVL[_k])[_distanceKey]; 815 } 816 817 return _distances; 818} 819 820