src/jalview/analysis/NJTree.java

   1 /*
   2  * Jalview - A Sequence Alignment Editor and Viewer ($$Version-Rel$$)
   3  * Copyright (C) $$Year-Rel$$ The Jalview Authors
   4  *
   5  * This file is part of Jalview.
   6  *
   7  * Jalview is free software: you can redistribute it and/or
   8  * modify it under the terms of the GNU General Public License
   9  * as published by the Free Software Foundation, either version 3
  10  * of the License, or (at your option) any later version.
  11  *
  12  * Jalview is distributed in the hope that it will be useful, but
  13  * WITHOUT ANY WARRANTY; without even the implied warranty
  14  * of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  15  * PURPOSE.  See the GNU General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU General Public License
  18  * along with Jalview.  If not, see <http://www.gnu.org/licenses/>.
  19  * The Jalview Authors are detailed in the 'AUTHORS' file.
  20  */
  21 package jalview.analysis;
  22
  23 import jalview.analysis.scoremodels.ScoreModels;
  24 import jalview.api.analysis.DistanceScoreModelI;
  25 import jalview.api.analysis.ScoreModelI;
  26 import jalview.api.analysis.SimilarityScoreModelI;
  27 import jalview.datamodel.AlignmentView;
  28 import jalview.datamodel.BinaryNode;
  29 import jalview.datamodel.CigarArray;
  30 import jalview.datamodel.NodeTransformI;
  31 import jalview.datamodel.SeqCigar;
  32 import jalview.datamodel.Sequence;
  33 import jalview.datamodel.SequenceI;
  34 import jalview.datamodel.SequenceNode;
  35 import jalview.io.NewickFile;
  36 import jalview.math.MatrixI;
  37
  38 import java.util.Enumeration;
  39 import java.util.List;
  40 import java.util.Vector;
  41
  42 /**
  43  * DOCUMENT ME!
  44  *
  45  * @author $author$
  46  * @version $Revision$
  47  */
  48 public class NJTree
  49 {
  50   /*
  51    * 'methods'
  52    */
  53   public static final String AVERAGE_DISTANCE = "AV";
  54
  55   public static final String NEIGHBOUR_JOINING = "NJ";
  56
  57   public static final String FROM_FILE = "FromFile";
  58
  59   Vector<Cluster> cluster;
  60
  61   SequenceI[] sequence;
  62
  63   // SequenceData is a string representation of what the user
  64   // sees. The display may contain hidden columns.
  65   public AlignmentView seqData = null;
  66
  67   int[] done;
  68
  69   int noseqs;
  70
  71   int noClus;
  72
  73   MatrixI distance;
  74
  75   int mini;
  76
  77   int minj;
  78
  79   double ri;
  80
  81   double rj;
  82
  83   Vector<SequenceNode> groups = new Vector<SequenceNode>();
  84
  85   SequenceNode maxdist;
  86
  87   SequenceNode top;
  88
  89   double maxDistValue;
  90
  91   double maxheight;
  92
  93   int ycount;
  94
  95   Vector<SequenceNode> node;
  96
  97   String type;
  98
  99   String pwtype;
 100
 101   Object found = null;
 102
 103   boolean hasDistances = true; // normal case for jalview trees
 104
 105   boolean hasBootstrap = false; // normal case for jalview trees
 106
 107   private boolean hasRootDistance = true;
 108
 109   /**
 110    * Create a new NJTree object with leaves associated with sequences in seqs,
 111    * and original alignment data represented by Cigar strings.
 112    *
 113    * @param seqs
 114    *          SequenceI[]
 115    * @param odata
 116    *          Cigar[]
 117    * @param treefile
 118    *          NewickFile
 119    */
 120   public NJTree(SequenceI[] seqs, AlignmentView odata, NewickFile treefile)
 121   {
 122     this(seqs, treefile);
 123     if (odata != null)
 124     {
 125       seqData = odata;
 126     }
 127     /*
 128      * sequenceString = new String[odata.length]; char gapChar =
 129      * jalview.util.Comparison.GapChars.charAt(0); for (int i = 0; i <
 130      * odata.length; i++) { SequenceI oseq_aligned = odata[i].getSeq(gapChar);
 131      * sequenceString[i] = oseq_aligned.getSequence(); }
 132      */
 133   }
 134
 135   /**
 136    * Creates a new NJTree object from a tree from an external source
 137    *
 138    * @param seqs
 139    *          SequenceI which should be associated with leafs of treefile
 140    * @param treefile
 141    *          A parsed tree
 142    */
 143   public NJTree(SequenceI[] seqs, NewickFile treefile)
 144   {
 145     this.sequence = seqs;
 146     top = treefile.getTree();
 147
 148     /**
 149      * There is no dependent alignment to be recovered from an imported tree.
 150      *
 151      * if (sequenceString == null) { sequenceString = new String[seqs.length];
 152      * for (int i = 0; i < seqs.length; i++) { sequenceString[i] =
 153      * seqs[i].getSequence(); } }
 154      */
 155
 156     hasDistances = treefile.HasDistances();
 157     hasBootstrap = treefile.HasBootstrap();
 158     hasRootDistance = treefile.HasRootDistance();
 159
 160     maxheight = findHeight(top);
 161
 162     SequenceIdMatcher algnIds = new SequenceIdMatcher(seqs);
 163
 164     Vector<SequenceNode> leaves = findLeaves(top);
 165
 166     int i = 0;
 167     int namesleft = seqs.length;
 168
 169     SequenceNode j;
 170     SequenceI nam;
 171     String realnam;
 172     Vector<SequenceI> one2many = new Vector<SequenceI>();
 173     int countOne2Many = 0;
 174     while (i < leaves.size())
 175     {
 176       j = leaves.elementAt(i++);
 177       realnam = j.getName();
 178       nam = null;
 179
 180       if (namesleft > -1)
 181       {
 182         nam = algnIds.findIdMatch(realnam);
 183       }
 184
 185       if (nam != null)
 186       {
 187         j.setElement(nam);
 188         if (one2many.contains(nam))
 189         {
 190           countOne2Many++;
 191           // if (jalview.bin.Cache.log.isDebugEnabled())
 192           // jalview.bin.Cache.log.debug("One 2 many relationship for
 193           // "+nam.getName());
 194         }
 195         else
 196         {
 197           one2many.addElement(nam);
 198           namesleft--;
 199         }
 200       }
 201       else
 202       {
 203         j.setElement(new Sequence(realnam, "THISISAPLACEHLDER"));
 204         j.setPlaceholder(true);
 205       }
 206     }
 207     // if (jalview.bin.Cache.log.isDebugEnabled() && countOne2Many>0) {
 208     // jalview.bin.Cache.log.debug("There were "+countOne2Many+" alignment
 209     // sequence ids (out of "+one2many.size()+" unique ids) linked to two or
 210     // more leaves.");
 211     // }
 212     // one2many.clear();
 213   }
 214
 215   /**
 216    * Creates a new NJTree object.
 217    *
 218    * @param sequence
 219    *          DOCUMENT ME!
 220    * @param treeType
 221    *          DOCUMENT ME!
 222    * @param modelType
 223    *          DOCUMENT ME!
 224    * @param start
 225    *          DOCUMENT ME!
 226    * @param end
 227    *          DOCUMENT ME!
 228    */
 229   public NJTree(SequenceI[] sqs, AlignmentView seqView, String treeType,
 230           String modelType, ScoreModelI sm, int start, int end)
 231   {
 232     this.sequence = sqs;
 233     this.node = new Vector<SequenceNode>();
 234     this.type = treeType;
 235     this.pwtype = modelType;
 236     if (seqView != null)
 237     {
 238       this.seqData = seqView;
 239     }
 240     else
 241     {
 242       SeqCigar[] seqs = new SeqCigar[sequence.length];
 243       for (int i = 0; i < sequence.length; i++)
 244       {
 245         seqs[i] = new SeqCigar(sequence[i], start, end);
 246       }
 247       CigarArray sdata = new CigarArray(seqs);
 248       sdata.addOperation(CigarArray.M, end - start + 1);
 249       this.seqData = new AlignmentView(sdata, start);
 250     }
 251     // System.err.println("Made seqData");// dbg
 252     if (!(treeType.equals(NEIGHBOUR_JOINING)))
 253     {
 254       treeType = AVERAGE_DISTANCE;
 255     }
 256
 257     if (sm == null && !(modelType.equals("PID")))
 258     {
 259       if (ScoreModels.getInstance().forName(modelType) == null)
 260       {
 261         modelType = "BLOSUM62";
 262       }
 263     }
 264
 265     int i = 0;
 266
 267     done = new int[sequence.length];
 268
 269     while ((i < sequence.length) && (sequence[i] != null))
 270     {
 271       done[i] = 0;
 272       i++;
 273     }
 274
 275     noseqs = i++;
 276
 277     if (sm instanceof DistanceScoreModelI)
 278     {
 279       distance = ((DistanceScoreModelI) sm).findDistances(seqData);
 280     }
 281     else if (sm instanceof SimilarityScoreModelI)
 282     {
 283       /*
 284        * compute similarity and invert it to give a distance measure
 285        */
 286       MatrixI result = ((SimilarityScoreModelI) sm)
 287               .findSimilarities(seqData);
 288       double maxScore = result.getMaxValue();
 289       result.subtractAllFrom(maxScore);
 290       distance = result;
 291     }
 292
 293     // System.err.println("Made distances");// dbg
 294     makeLeaves();
 295     // System.err.println("Made leaves");// dbg
 296
 297     noClus = cluster.size();
 298
 299     cluster();
 300     // System.err.println("Made clusters");// dbg
 301
 302   }
 303
 304   /**
 305    * Generate a string representation of the Tree
 306    *
 307    * @return Newick File with all tree data available
 308    */
 309   @Override
 310   public String toString()
 311   {
 312     jalview.io.NewickFile fout = new jalview.io.NewickFile(getTopNode());
 313
 314     return fout.print(isHasBootstrap(), isHasDistances(),
 315             isHasRootDistance()); // output all data available for tree
 316   }
 317
 318   /**
 319    *
 320    * used when the alignment associated to a tree has changed.
 321    *
 322    * @param list
 323    *          Sequence set to be associated with tree nodes
 324    */
 325   public void UpdatePlaceHolders(List<SequenceI> list)
 326   {
 327     Vector<SequenceNode> leaves = findLeaves(top);
 328
 329     int sz = leaves.size();
 330     SequenceIdMatcher seqmatcher = null;
 331     int i = 0;
 332
 333     while (i < sz)
 334     {
 335       SequenceNode leaf = leaves.elementAt(i++);
 336
 337       if (list.contains(leaf.element()))
 338       {
 339         leaf.setPlaceholder(false);
 340       }
 341       else
 342       {
 343         if (seqmatcher == null)
 344         {
 345           // Only create this the first time we need it
 346           SequenceI[] seqs = new SequenceI[list.size()];
 347
 348           for (int j = 0; j < seqs.length; j++)
 349           {
 350             seqs[j] = list.get(j);
 351           }
 352
 353           seqmatcher = new SequenceIdMatcher(seqs);
 354         }
 355
 356         SequenceI nam = seqmatcher.findIdMatch(leaf.getName());
 357
 358         if (nam != null)
 359         {
 360           if (!leaf.isPlaceholder())
 361           {
 362             // remapping the node to a new sequenceI - should remove any refs to
 363             // old one.
 364             // TODO - make many sequenceI to one leaf mappings possible!
 365             // (JBPNote)
 366           }
 367           leaf.setPlaceholder(false);
 368           leaf.setElement(nam);
 369         }
 370         else
 371         {
 372           if (!leaf.isPlaceholder())
 373           {
 374             // Construct a new placeholder sequence object for this leaf
 375             leaf.setElement(new Sequence(leaf.getName(),
 376                     "THISISAPLACEHLDER"));
 377           }
 378           leaf.setPlaceholder(true);
 379
 380         }
 381       }
 382     }
 383   }
 384
 385   /**
 386    * rename any nodes according to their associated sequence. This will modify
 387    * the tree's metadata! (ie the original NewickFile or newly generated
 388    * BinaryTree's label data)
 389    */
 390   public void renameAssociatedNodes()
 391   {
 392     applyToNodes(new NodeTransformI()
 393     {
 394
 395       @Override
 396       public void transform(BinaryNode nd)
 397       {
 398         Object el = nd.element();
 399         if (el != null && el instanceof SequenceI)
 400         {
 401           nd.setName(((SequenceI) el).getName());
 402         }
 403       }
 404     });
 405   }
 406
 407   /**
 408    * DOCUMENT ME!
 409    */
 410   public void cluster()
 411   {
 412     while (noClus > 2)
 413     {
 414       if (type.equals(NEIGHBOUR_JOINING))
 415       {
 416         findMinNJDistance();
 417       }
 418       else
 419       {
 420         findMinDistance();
 421       }
 422
 423       Cluster c = joinClusters(mini, minj);
 424
 425       done[minj] = 1;
 426
 427       cluster.setElementAt(null, minj);
 428       cluster.setElementAt(c, mini);
 429
 430       noClus--;
 431     }
 432
 433     boolean onefound = false;
 434
 435     int one = -1;
 436     int two = -1;
 437
 438     for (int i = 0; i < noseqs; i++)
 439     {
 440       if (done[i] != 1)
 441       {
 442         if (onefound == false)
 443         {
 444           two = i;
 445           onefound = true;
 446         }
 447         else
 448         {
 449           one = i;
 450         }
 451       }
 452     }
 453
 454     joinClusters(one, two);
 455     top = (node.elementAt(one));
 456
 457     reCount(top);
 458     findHeight(top);
 459     findMaxDist(top);
 460   }
 461
 462   /**
 463    * DOCUMENT ME!
 464    *
 465    * @param i
 466    *          DOCUMENT ME!
 467    * @param j
 468    *          DOCUMENT ME!
 469    *
 470    * @return DOCUMENT ME!
 471    */
 472   public Cluster joinClusters(int i, int j)
 473   {
 474     double dist = distance.getValue(i, j);
 475
 476     int noi = cluster.elementAt(i).value.length;
 477     int noj = cluster.elementAt(j).value.length;
 478
 479     int[] value = new int[noi + noj];
 480
 481     for (int ii = 0; ii < noi; ii++)
 482     {
 483       value[ii] = cluster.elementAt(i).value[ii];
 484     }
 485
 486     for (int ii = noi; ii < (noi + noj); ii++)
 487     {
 488       value[ii] = cluster.elementAt(j).value[ii - noi];
 489     }
 490
 491     Cluster c = new Cluster(value);
 492
 493     ri = findr(i, j);
 494     rj = findr(j, i);
 495
 496     if (type.equals(NEIGHBOUR_JOINING))
 497     {
 498       findClusterNJDistance(i, j);
 499     }
 500     else
 501     {
 502       findClusterDistance(i, j);
 503     }
 504
 505     SequenceNode sn = new SequenceNode();
 506
 507     sn.setLeft((node.elementAt(i)));
 508     sn.setRight((node.elementAt(j)));
 509
 510     SequenceNode tmpi = (node.elementAt(i));
 511     SequenceNode tmpj = (node.elementAt(j));
 512
 513     if (type.equals(NEIGHBOUR_JOINING))
 514     {
 515       findNewNJDistances(tmpi, tmpj, dist);
 516     }
 517     else
 518     {
 519       findNewDistances(tmpi, tmpj, dist);
 520     }
 521
 522     tmpi.setParent(sn);
 523     tmpj.setParent(sn);
 524
 525     node.setElementAt(sn, i);
 526
 527     return c;
 528   }
 529
 530   /**
 531    * DOCUMENT ME!
 532    *
 533    * @param tmpi
 534    *          DOCUMENT ME!
 535    * @param tmpj
 536    *          DOCUMENT ME!
 537    * @param dist
 538    *          DOCUMENT ME!
 539    */
 540   public void findNewNJDistances(SequenceNode tmpi, SequenceNode tmpj,
 541           double dist)
 542   {
 543
 544     tmpi.dist = ((dist + ri) - rj) / 2;
 545     tmpj.dist = (dist - tmpi.dist);
 546
 547     if (tmpi.dist < 0)
 548     {
 549       tmpi.dist = 0;
 550     }
 551
 552     if (tmpj.dist < 0)
 553     {
 554       tmpj.dist = 0;
 555     }
 556   }
 557
 558   /**
 559    * DOCUMENT ME!
 560    *
 561    * @param tmpi
 562    *          DOCUMENT ME!
 563    * @param tmpj
 564    *          DOCUMENT ME!
 565    * @param dist
 566    *          DOCUMENT ME!
 567    */
 568   public void findNewDistances(SequenceNode tmpi, SequenceNode tmpj,
 569           double dist)
 570   {
 571     double ih = 0;
 572     double jh = 0;
 573
 574     SequenceNode sni = tmpi;
 575     SequenceNode snj = tmpj;
 576
 577     while (sni != null)
 578     {
 579       ih = ih + sni.dist;
 580       sni = (SequenceNode) sni.left();
 581     }
 582
 583     while (snj != null)
 584     {
 585       jh = jh + snj.dist;
 586       snj = (SequenceNode) snj.left();
 587     }
 588
 589     tmpi.dist = ((dist / 2) - ih);
 590     tmpj.dist = ((dist / 2) - jh);
 591   }
 592
 593   /**
 594    * DOCUMENT ME!
 595    *
 596    * @param i
 597    *          DOCUMENT ME!
 598    * @param j
 599    *          DOCUMENT ME!
 600    */
 601   public void findClusterDistance(int i, int j)
 602   {
 603     int noi = cluster.elementAt(i).value.length;
 604     int noj = cluster.elementAt(j).value.length;
 605
 606     // New distances from cluster to others
 607     double[] newdist = new double[noseqs];
 608
 609     for (int l = 0; l < noseqs; l++)
 610     {
 611       if ((l != i) && (l != j))
 612       {
 613         // newdist[l] = ((distance[i][l] * noi) + (distance[j][l] * noj))
 614         // / (noi + noj);
 615         newdist[l] = ((distance.getValue(i, l) * noi) + (distance.getValue(
 616                 j, l) * noj))
 617                 / (noi + noj);
 618       }
 619       else
 620       {
 621         newdist[l] = 0;
 622       }
 623     }
 624
 625     for (int ii = 0; ii < noseqs; ii++)
 626     {
 627       // distance[i][ii] = newdist[ii];
 628       // distance[ii][i] = newdist[ii];
 629       distance.setValue(i, ii, newdist[ii]);
 630       distance.setValue(ii, i, newdist[ii]);
 631     }
 632   }
 633
 634   /**
 635    * DOCUMENT ME!
 636    *
 637    * @param i
 638    *          DOCUMENT ME!
 639    * @param j
 640    *          DOCUMENT ME!
 641    */
 642   public void findClusterNJDistance(int i, int j)
 643   {
 644
 645     // New distances from cluster to others
 646     double[] newdist = new double[noseqs];
 647
 648     for (int l = 0; l < noseqs; l++)
 649     {
 650       if ((l != i) && (l != j))
 651       {
 652         // newdist[l] = ((distance[i][l] + distance[j][l]) - distance[i][j]) /
 653         // 2;
 654         newdist[l] = (distance.getValue(i, l) + distance.getValue(j, l) - distance
 655                 .getValue(i, j)) / 2;
 656       }
 657       else
 658       {
 659         newdist[l] = 0;
 660       }
 661     }
 662
 663     for (int ii = 0; ii < noseqs; ii++)
 664     {
 665       // distance[i][ii] = newdist[ii];
 666       // distance[ii][i] = newdist[ii];
 667       distance.setValue(i, ii, newdist[ii]);
 668       distance.setValue(ii, i, newdist[ii]);
 669     }
 670   }
 671
 672   /**
 673    * DOCUMENT ME!
 674    *
 675    * @param i
 676    *          DOCUMENT ME!
 677    * @param j
 678    *          DOCUMENT ME!
 679    *
 680    * @return DOCUMENT ME!
 681    */
 682   public double findr(int i, int j)
 683   {
 684     double tmp = 1;
 685
 686     for (int k = 0; k < noseqs; k++)
 687     {
 688       if ((k != i) && (k != j) && (done[k] != 1))
 689       {
 690         // tmp = tmp + distance[i][k];
 691         tmp = tmp + distance.getValue(i, k);
 692       }
 693     }
 694
 695     if (noClus > 2)
 696     {
 697       tmp = tmp / (noClus - 2);
 698     }
 699
 700     return tmp;
 701   }
 702
 703   /**
 704    * DOCUMENT ME!
 705    *
 706    * @return DOCUMENT ME!
 707    */
 708   public double findMinNJDistance()
 709   {
 710     double min = Double.MAX_VALUE;
 711
 712     for (int i = 0; i < (noseqs - 1); i++)
 713     {
 714       for (int j = i + 1; j < noseqs; j++)
 715       {
 716         if ((done[i] != 1) && (done[j] != 1))
 717         {
 718           // float tmp = distance[i][j] - (findr(i, j) + findr(j, i));
 719           double tmp = distance.getValue(i, j)
 720                   - (findr(i, j) + findr(j, i));
 721
 722           if (tmp < min)
 723           {
 724             mini = i;
 725             minj = j;
 726
 727             min = tmp;
 728           }
 729         }
 730       }
 731     }
 732
 733     return min;
 734   }
 735
 736   /**
 737    * DOCUMENT ME!
 738    *
 739    * @return DOCUMENT ME!
 740    */
 741   public double findMinDistance()
 742   {
 743     double min = Double.MAX_VALUE;
 744
 745     for (int i = 0; i < (noseqs - 1); i++)
 746     {
 747       for (int j = i + 1; j < noseqs; j++)
 748       {
 749         if ((done[i] != 1) && (done[j] != 1))
 750         {
 751           // if (distance[i][j] < min)
 752           if (distance.getValue(i, j) < min)
 753           {
 754             mini = i;
 755             minj = j;
 756
 757             // min = distance[i][j];
 758             min = distance.getValue(i, j);
 759           }
 760         }
 761       }
 762     }
 763
 764     return min;
 765   }
 766
 767   /**
 768    * Calculate a distance matrix given the sequence input data and score model
 769    *
 770    * @return
 771    */
 772   public MatrixI findDistances(ScoreModelI scoreModel)
 773   {
 774     MatrixI result = null;
 775
 776     if (scoreModel == null)
 777     {
 778       // Resolve substitution model
 779       scoreModel = ScoreModels.getInstance().forName(pwtype);
 780       if (scoreModel == null)
 781       {
 782         scoreModel = ScoreModels.getInstance().forName("BLOSUM62");
 783       }
 784     }
 785     if (scoreModel instanceof DistanceScoreModelI)
 786     {
 787       result = ((DistanceScoreModelI) scoreModel).findDistances(seqData);
 788     }
 789     else if (scoreModel instanceof SimilarityScoreModelI)
 790     {
 791       /*
 792        * compute similarity and invert it to give a distance measure
 793        */
 794       result = ((SimilarityScoreModelI) scoreModel)
 795               .findSimilarities(seqData);
 796       double maxScore = result.getMaxValue();
 797       result.subtractAllFrom(maxScore);
 798     }
 799     else
 800     {
 801       System.err
 802               .println("Unexpected type of score model, can't compute distances");
 803     }
 804
 805     return result;
 806   }
 807
 808   /**
 809    * DOCUMENT ME!
 810    */
 811   public void makeLeaves()
 812   {
 813     cluster = new Vector<Cluster>();
 814
 815     for (int i = 0; i < noseqs; i++)
 816     {
 817       SequenceNode sn = new SequenceNode();
 818
 819       sn.setElement(sequence[i]);
 820       sn.setName(sequence[i].getName());
 821       node.addElement(sn);
 822
 823       int[] value = new int[1];
 824       value[0] = i;
 825
 826       Cluster c = new Cluster(value);
 827       cluster.addElement(c);
 828     }
 829   }
 830
 831   /**
 832    * Search for leaf nodes below (or at) the given node
 833    *
 834    * @param nd
 835    *          root node to search from
 836    *
 837    * @return
 838    */
 839   public Vector<SequenceNode> findLeaves(SequenceNode nd)
 840   {
 841     Vector<SequenceNode> leaves = new Vector<SequenceNode>();
 842     findLeaves(nd, leaves);
 843     return leaves;
 844   }
 845
 846   /**
 847    * Search for leaf nodes.
 848    *
 849    * @param nd
 850    *          root node to search from
 851    * @param leaves
 852    *          Vector of leaves to add leaf node objects too.
 853    *
 854    * @return Vector of leaf nodes on binary tree
 855    */
 856   Vector<SequenceNode> findLeaves(SequenceNode nd,
 857           Vector<SequenceNode> leaves)
 858   {
 859     if (nd == null)
 860     {
 861       return leaves;
 862     }
 863
 864     if ((nd.left() == null) && (nd.right() == null)) // Interior node
 865     // detection
 866     {
 867       leaves.addElement(nd);
 868
 869       return leaves;
 870     }
 871     else
 872     {
 873       /*
 874        * TODO: Identify internal nodes... if (node.isSequenceLabel()) {
 875        * leaves.addElement(node); }
 876        */
 877       findLeaves((SequenceNode) nd.left(), leaves);
 878       findLeaves((SequenceNode) nd.right(), leaves);
 879     }
 880
 881     return leaves;
 882   }
 883
 884   /**
 885    * Find the leaf node with a particular ycount
 886    *
 887    * @param nd
 888    *          initial point on tree to search from
 889    * @param count
 890    *          value to search for
 891    *
 892    * @return null or the node with ycound=count
 893    */
 894   public Object findLeaf(SequenceNode nd, int count)
 895   {
 896     found = _findLeaf(nd, count);
 897
 898     return found;
 899   }
 900
 901   /*
 902    * #see findLeaf(SequenceNode node, count)
 903    */
 904   public Object _findLeaf(SequenceNode nd, int count)
 905   {
 906     if (nd == null)
 907     {
 908       return null;
 909     }
 910
 911     if (nd.ycount == count)
 912     {
 913       found = nd.element();
 914
 915       return found;
 916     }
 917     else
 918     {
 919       _findLeaf((SequenceNode) nd.left(), count);
 920       _findLeaf((SequenceNode) nd.right(), count);
 921     }
 922
 923     return found;
 924   }
 925
 926   /**
 927    * printNode is mainly for debugging purposes.
 928    *
 929    * @param nd
 930    *          SequenceNode
 931    */
 932   public void printNode(SequenceNode nd)
 933   {
 934     if (nd == null)
 935     {
 936       return;
 937     }
 938
 939     if ((nd.left() == null) && (nd.right() == null))
 940     {
 941       System.out.println("Leaf = " + ((SequenceI) nd.element()).getName());
 942       System.out.println("Dist " + nd.dist);
 943       System.out.println("Boot " + nd.getBootstrap());
 944     }
 945     else
 946     {
 947       System.out.println("Dist " + nd.dist);
 948       printNode((SequenceNode) nd.left());
 949       printNode((SequenceNode) nd.right());
 950     }
 951   }
 952
 953   /**
 954    * DOCUMENT ME!
 955    *
 956    * @param nd
 957    *          DOCUMENT ME!
 958    */
 959   public void findMaxDist(SequenceNode nd)
 960   {
 961     if (nd == null)
 962     {
 963       return;
 964     }
 965
 966     if ((nd.left() == null) && (nd.right() == null))
 967     {
 968       double dist = nd.dist;
 969
 970       if (dist > maxDistValue)
 971       {
 972         maxdist = nd;
 973         maxDistValue = dist;
 974       }
 975     }
 976     else
 977     {
 978       findMaxDist((SequenceNode) nd.left());
 979       findMaxDist((SequenceNode) nd.right());
 980     }
 981   }
 982
 983   /**
 984    * DOCUMENT ME!
 985    *
 986    * @return DOCUMENT ME!
 987    */
 988   public Vector<SequenceNode> getGroups()
 989   {
 990     return groups;
 991   }
 992
 993   /**
 994    * DOCUMENT ME!
 995    *
 996    * @return DOCUMENT ME!
 997    */
 998   public double getMaxHeight()
 999   {
1000     return maxheight;
1001   }
1002
1003   /**
1004    * DOCUMENT ME!
1005    *
1006    * @param nd
1007    *          DOCUMENT ME!
1008    * @param threshold
1009    *          DOCUMENT ME!
1010    */
1011   public void groupNodes(SequenceNode nd, float threshold)
1012   {
1013     if (nd == null)
1014     {
1015       return;
1016     }
1017
1018     if ((nd.height / maxheight) > threshold)
1019     {
1020       groups.addElement(nd);
1021     }
1022     else
1023     {
1024       groupNodes((SequenceNode) nd.left(), threshold);
1025       groupNodes((SequenceNode) nd.right(), threshold);
1026     }
1027   }
1028
1029   /**
1030    * DOCUMENT ME!
1031    *
1032    * @param nd
1033    *          DOCUMENT ME!
1034    *
1035    * @return DOCUMENT ME!
1036    */
1037   public double findHeight(SequenceNode nd)
1038   {
1039     if (nd == null)
1040     {
1041       return maxheight;
1042     }
1043
1044     if ((nd.left() == null) && (nd.right() == null))
1045     {
1046       nd.height = ((SequenceNode) nd.parent()).height + nd.dist;
1047
1048       if (nd.height > maxheight)
1049       {
1050         return nd.height;
1051       }
1052       else
1053       {
1054         return maxheight;
1055       }
1056     }
1057     else
1058     {
1059       if (nd.parent() != null)
1060       {
1061         nd.height = ((SequenceNode) nd.parent()).height + nd.dist;
1062       }
1063       else
1064       {
1065         maxheight = 0;
1066         nd.height = (float) 0.0;
1067       }
1068
1069       maxheight = findHeight((SequenceNode) (nd.left()));
1070       maxheight = findHeight((SequenceNode) (nd.right()));
1071     }
1072
1073     return maxheight;
1074   }
1075
1076   /**
1077    * DOCUMENT ME!
1078    *
1079    * @return DOCUMENT ME!
1080    */
1081   public SequenceNode reRoot()
1082   {
1083     if (maxdist != null)
1084     {
1085       ycount = 0;
1086
1087       double tmpdist = maxdist.dist;
1088
1089       // New top
1090       SequenceNode sn = new SequenceNode();
1091       sn.setParent(null);
1092
1093       // New right hand of top
1094       SequenceNode snr = (SequenceNode) maxdist.parent();
1095       changeDirection(snr, maxdist);
1096       System.out.println("Printing reversed tree");
1097       printN(snr);
1098       snr.dist = tmpdist / 2;
1099       maxdist.dist = tmpdist / 2;
1100
1101       snr.setParent(sn);
1102       maxdist.setParent(sn);
1103
1104       sn.setRight(snr);
1105       sn.setLeft(maxdist);
1106
1107       top = sn;
1108
1109       ycount = 0;
1110       reCount(top);
1111       findHeight(top);
1112     }
1113
1114     return top;
1115   }
1116
1117   /**
1118    *
1119    * @return true if original sequence data can be recovered
1120    */
1121   public boolean hasOriginalSequenceData()
1122   {
1123     return seqData != null;
1124   }
1125
1126   /**
1127    * Returns original alignment data used for calculation - or null where not
1128    * available.
1129    *
1130    * @return null or cut'n'pasteable alignment
1131    */
1132   public String printOriginalSequenceData(char gapChar)
1133   {
1134     if (seqData == null)
1135     {
1136       return null;
1137     }
1138
1139     StringBuffer sb = new StringBuffer();
1140     String[] seqdatas = seqData.getSequenceStrings(gapChar);
1141     for (int i = 0; i < seqdatas.length; i++)
1142     {
1143       sb.append(new jalview.util.Format("%-" + 15 + "s").form(sequence[i]
1144               .getName()));
1145       sb.append(" " + seqdatas[i] + "\n");
1146     }
1147     return sb.toString();
1148   }
1149
1150   /**
1151    * DOCUMENT ME!
1152    *
1153    * @param nd
1154    *          DOCUMENT ME!
1155    */
1156   public void printN(SequenceNode nd)
1157   {
1158     if (nd == null)
1159     {
1160       return;
1161     }
1162
1163     if ((nd.left() != null) && (nd.right() != null))
1164     {
1165       printN((SequenceNode) nd.left());
1166       printN((SequenceNode) nd.right());
1167     }
1168     else
1169     {
1170       System.out.println(" name = " + ((SequenceI) nd.element()).getName());
1171     }
1172
1173     System.out.println(" dist = " + nd.dist + " " + nd.count + " "
1174             + nd.height);
1175   }
1176
1177   /**
1178    * DOCUMENT ME!
1179    *
1180    * @param nd
1181    *          DOCUMENT ME!
1182    */
1183   public void reCount(SequenceNode nd)
1184   {
1185     ycount = 0;
1186     _lycount = 0;
1187     // _lylimit = this.node.size();
1188     _reCount(nd);
1189   }
1190
1191   private long _lycount = 0, _lylimit = 0;
1192
1193   /**
1194    * DOCUMENT ME!
1195    *
1196    * @param nd
1197    *          DOCUMENT ME!
1198    */
1199   public void _reCount(SequenceNode nd)
1200   {
1201     // if (_lycount<_lylimit)
1202     // {
1203     // System.err.println("Warning: depth of _recount greater than number of nodes.");
1204     // }
1205     if (nd == null)
1206     {
1207       return;
1208     }
1209     _lycount++;
1210
1211     if ((nd.left() != null) && (nd.right() != null))
1212     {
1213
1214       _reCount((SequenceNode) nd.left());
1215       _reCount((SequenceNode) nd.right());
1216
1217       SequenceNode l = (SequenceNode) nd.left();
1218       SequenceNode r = (SequenceNode) nd.right();
1219
1220       nd.count = l.count + r.count;
1221       nd.ycount = (l.ycount + r.ycount) / 2;
1222     }
1223     else
1224     {
1225       nd.count = 1;
1226       nd.ycount = ycount++;
1227     }
1228     _lycount--;
1229   }
1230
1231   /**
1232    * DOCUMENT ME!
1233    *
1234    * @param nd
1235    *          DOCUMENT ME!
1236    */
1237   public void swapNodes(SequenceNode nd)
1238   {
1239     if (nd == null)
1240     {
1241       return;
1242     }
1243
1244     SequenceNode tmp = (SequenceNode) nd.left();
1245
1246     nd.setLeft(nd.right());
1247     nd.setRight(tmp);
1248   }
1249
1250   /**
1251    * DOCUMENT ME!
1252    *
1253    * @param nd
1254    *          DOCUMENT ME!
1255    * @param dir
1256    *          DOCUMENT ME!
1257    */
1258   public void changeDirection(SequenceNode nd, SequenceNode dir)
1259   {
1260     if (nd == null)
1261     {
1262       return;
1263     }
1264
1265     if (nd.parent() != top)
1266     {
1267       changeDirection((SequenceNode) nd.parent(), nd);
1268
1269       SequenceNode tmp = (SequenceNode) nd.parent();
1270
1271       if (dir == nd.left())
1272       {
1273         nd.setParent(dir);
1274         nd.setLeft(tmp);
1275       }
1276       else if (dir == nd.right())
1277       {
1278         nd.setParent(dir);
1279         nd.setRight(tmp);
1280       }
1281     }
1282     else
1283     {
1284       if (dir == nd.left())
1285       {
1286         nd.setParent(nd.left());
1287
1288         if (top.left() == nd)
1289         {
1290           nd.setRight(top.right());
1291         }
1292         else
1293         {
1294           nd.setRight(top.left());
1295         }
1296       }
1297       else
1298       {
1299         nd.setParent(nd.right());
1300
1301         if (top.left() == nd)
1302         {
1303           nd.setLeft(top.right());
1304         }
1305         else
1306         {
1307           nd.setLeft(top.left());
1308         }
1309       }
1310     }
1311   }
1312
1313   /**
1314    * DOCUMENT ME!
1315    *
1316    * @return DOCUMENT ME!
1317    */
1318   public SequenceNode getMaxDist()
1319   {
1320     return maxdist;
1321   }
1322
1323   /**
1324    * DOCUMENT ME!
1325    *
1326    * @return DOCUMENT ME!
1327    */
1328   public SequenceNode getTopNode()
1329   {
1330     return top;
1331   }
1332
1333   /**
1334    *
1335    * @return true if tree has real distances
1336    */
1337   public boolean isHasDistances()
1338   {
1339     return hasDistances;
1340   }
1341
1342   /**
1343    *
1344    * @return true if tree has real bootstrap values
1345    */
1346   public boolean isHasBootstrap()
1347   {
1348     return hasBootstrap;
1349   }
1350
1351   public boolean isHasRootDistance()
1352   {
1353     return hasRootDistance;
1354   }
1355
1356   /**
1357    * apply the given transform to all the nodes in the tree.
1358    *
1359    * @param nodeTransformI
1360    */
1361   public void applyToNodes(NodeTransformI nodeTransformI)
1362   {
1363     for (Enumeration<SequenceNode> nodes = node.elements(); nodes
1364             .hasMoreElements(); nodeTransformI.transform(nodes
1365             .nextElement()))
1366     {
1367       ;
1368     }
1369   }
1370 }
1371
1372 /**
1373  * DOCUMENT ME!
1374  *
1375  * @author $author$
1376  * @version $Revision$
1377  */
1378 // TODO what does this class have that int[] doesn't have already?
1379 class Cluster
1380 {
1381   int[] value;
1382
1383   /**
1384    * Creates a new Cluster object.
1385    *
1386    * @param value
1387    *          DOCUMENT ME!
1388    */
1389   public Cluster(int[] value)
1390   {
1391     this.value = value;
1392   }
1393 }