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       result.reverseRange(true);
 289       distance = result;
 290     }
 291
 292     // System.err.println("Made distances");// dbg
 293     makeLeaves();
 294     // System.err.println("Made leaves");// dbg
 295
 296     noClus = cluster.size();
 297
 298     cluster();
 299     // System.err.println("Made clusters");// dbg
 300
 301   }
 302
 303   /**
 304    * Generate a string representation of the Tree
 305    *
 306    * @return Newick File with all tree data available
 307    */
 308   @Override
 309   public String toString()
 310   {
 311     jalview.io.NewickFile fout = new jalview.io.NewickFile(getTopNode());
 312
 313     return fout.print(isHasBootstrap(), isHasDistances(),
 314             isHasRootDistance()); // output all data available for tree
 315   }
 316
 317   /**
 318    *
 319    * used when the alignment associated to a tree has changed.
 320    *
 321    * @param list
 322    *          Sequence set to be associated with tree nodes
 323    */
 324   public void UpdatePlaceHolders(List<SequenceI> list)
 325   {
 326     Vector<SequenceNode> leaves = findLeaves(top);
 327
 328     int sz = leaves.size();
 329     SequenceIdMatcher seqmatcher = null;
 330     int i = 0;
 331
 332     while (i < sz)
 333     {
 334       SequenceNode leaf = leaves.elementAt(i++);
 335
 336       if (list.contains(leaf.element()))
 337       {
 338         leaf.setPlaceholder(false);
 339       }
 340       else
 341       {
 342         if (seqmatcher == null)
 343         {
 344           // Only create this the first time we need it
 345           SequenceI[] seqs = new SequenceI[list.size()];
 346
 347           for (int j = 0; j < seqs.length; j++)
 348           {
 349             seqs[j] = list.get(j);
 350           }
 351
 352           seqmatcher = new SequenceIdMatcher(seqs);
 353         }
 354
 355         SequenceI nam = seqmatcher.findIdMatch(leaf.getName());
 356
 357         if (nam != null)
 358         {
 359           if (!leaf.isPlaceholder())
 360           {
 361             // remapping the node to a new sequenceI - should remove any refs to
 362             // old one.
 363             // TODO - make many sequenceI to one leaf mappings possible!
 364             // (JBPNote)
 365           }
 366           leaf.setPlaceholder(false);
 367           leaf.setElement(nam);
 368         }
 369         else
 370         {
 371           if (!leaf.isPlaceholder())
 372           {
 373             // Construct a new placeholder sequence object for this leaf
 374             leaf.setElement(new Sequence(leaf.getName(),
 375                     "THISISAPLACEHLDER"));
 376           }
 377           leaf.setPlaceholder(true);
 378
 379         }
 380       }
 381     }
 382   }
 383
 384   /**
 385    * rename any nodes according to their associated sequence. This will modify
 386    * the tree's metadata! (ie the original NewickFile or newly generated
 387    * BinaryTree's label data)
 388    */
 389   public void renameAssociatedNodes()
 390   {
 391     applyToNodes(new NodeTransformI()
 392     {
 393
 394       @Override
 395       public void transform(BinaryNode nd)
 396       {
 397         Object el = nd.element();
 398         if (el != null && el instanceof SequenceI)
 399         {
 400           nd.setName(((SequenceI) el).getName());
 401         }
 402       }
 403     });
 404   }
 405
 406   /**
 407    * DOCUMENT ME!
 408    */
 409   public void cluster()
 410   {
 411     while (noClus > 2)
 412     {
 413       if (type.equals(NEIGHBOUR_JOINING))
 414       {
 415         findMinNJDistance();
 416       }
 417       else
 418       {
 419         findMinDistance();
 420       }
 421
 422       Cluster c = joinClusters(mini, minj);
 423
 424       done[minj] = 1;
 425
 426       cluster.setElementAt(null, minj);
 427       cluster.setElementAt(c, mini);
 428
 429       noClus--;
 430     }
 431
 432     boolean onefound = false;
 433
 434     int one = -1;
 435     int two = -1;
 436
 437     for (int i = 0; i < noseqs; i++)
 438     {
 439       if (done[i] != 1)
 440       {
 441         if (onefound == false)
 442         {
 443           two = i;
 444           onefound = true;
 445         }
 446         else
 447         {
 448           one = i;
 449         }
 450       }
 451     }
 452
 453     joinClusters(one, two);
 454     top = (node.elementAt(one));
 455
 456     reCount(top);
 457     findHeight(top);
 458     findMaxDist(top);
 459   }
 460
 461   /**
 462    * DOCUMENT ME!
 463    *
 464    * @param i
 465    *          DOCUMENT ME!
 466    * @param j
 467    *          DOCUMENT ME!
 468    *
 469    * @return DOCUMENT ME!
 470    */
 471   public Cluster joinClusters(int i, int j)
 472   {
 473     double dist = distance.getValue(i, j);
 474
 475     int noi = cluster.elementAt(i).value.length;
 476     int noj = cluster.elementAt(j).value.length;
 477
 478     int[] value = new int[noi + noj];
 479
 480     for (int ii = 0; ii < noi; ii++)
 481     {
 482       value[ii] = cluster.elementAt(i).value[ii];
 483     }
 484
 485     for (int ii = noi; ii < (noi + noj); ii++)
 486     {
 487       value[ii] = cluster.elementAt(j).value[ii - noi];
 488     }
 489
 490     Cluster c = new Cluster(value);
 491
 492     ri = findr(i, j);
 493     rj = findr(j, i);
 494
 495     if (type.equals(NEIGHBOUR_JOINING))
 496     {
 497       findClusterNJDistance(i, j);
 498     }
 499     else
 500     {
 501       findClusterDistance(i, j);
 502     }
 503
 504     SequenceNode sn = new SequenceNode();
 505
 506     sn.setLeft((node.elementAt(i)));
 507     sn.setRight((node.elementAt(j)));
 508
 509     SequenceNode tmpi = (node.elementAt(i));
 510     SequenceNode tmpj = (node.elementAt(j));
 511
 512     if (type.equals(NEIGHBOUR_JOINING))
 513     {
 514       findNewNJDistances(tmpi, tmpj, dist);
 515     }
 516     else
 517     {
 518       findNewDistances(tmpi, tmpj, dist);
 519     }
 520
 521     tmpi.setParent(sn);
 522     tmpj.setParent(sn);
 523
 524     node.setElementAt(sn, i);
 525
 526     return c;
 527   }
 528
 529   /**
 530    * DOCUMENT ME!
 531    *
 532    * @param tmpi
 533    *          DOCUMENT ME!
 534    * @param tmpj
 535    *          DOCUMENT ME!
 536    * @param dist
 537    *          DOCUMENT ME!
 538    */
 539   public void findNewNJDistances(SequenceNode tmpi, SequenceNode tmpj,
 540           double dist)
 541   {
 542
 543     tmpi.dist = ((dist + ri) - rj) / 2;
 544     tmpj.dist = (dist - tmpi.dist);
 545
 546     if (tmpi.dist < 0)
 547     {
 548       tmpi.dist = 0;
 549     }
 550
 551     if (tmpj.dist < 0)
 552     {
 553       tmpj.dist = 0;
 554     }
 555   }
 556
 557   /**
 558    * DOCUMENT ME!
 559    *
 560    * @param tmpi
 561    *          DOCUMENT ME!
 562    * @param tmpj
 563    *          DOCUMENT ME!
 564    * @param dist
 565    *          DOCUMENT ME!
 566    */
 567   public void findNewDistances(SequenceNode tmpi, SequenceNode tmpj,
 568           double dist)
 569   {
 570     double ih = 0;
 571     double jh = 0;
 572
 573     SequenceNode sni = tmpi;
 574     SequenceNode snj = tmpj;
 575
 576     while (sni != null)
 577     {
 578       ih = ih + sni.dist;
 579       sni = (SequenceNode) sni.left();
 580     }
 581
 582     while (snj != null)
 583     {
 584       jh = jh + snj.dist;
 585       snj = (SequenceNode) snj.left();
 586     }
 587
 588     tmpi.dist = ((dist / 2) - ih);
 589     tmpj.dist = ((dist / 2) - jh);
 590   }
 591
 592   /**
 593    * DOCUMENT ME!
 594    *
 595    * @param i
 596    *          DOCUMENT ME!
 597    * @param j
 598    *          DOCUMENT ME!
 599    */
 600   public void findClusterDistance(int i, int j)
 601   {
 602     int noi = cluster.elementAt(i).value.length;
 603     int noj = cluster.elementAt(j).value.length;
 604
 605     // New distances from cluster to others
 606     double[] newdist = new double[noseqs];
 607
 608     for (int l = 0; l < noseqs; l++)
 609     {
 610       if ((l != i) && (l != j))
 611       {
 612         // newdist[l] = ((distance[i][l] * noi) + (distance[j][l] * noj))
 613         // / (noi + noj);
 614         newdist[l] = ((distance.getValue(i, l) * noi) + (distance.getValue(
 615                 j, l) * noj))
 616                 / (noi + noj);
 617       }
 618       else
 619       {
 620         newdist[l] = 0;
 621       }
 622     }
 623
 624     for (int ii = 0; ii < noseqs; ii++)
 625     {
 626       // distance[i][ii] = newdist[ii];
 627       // distance[ii][i] = newdist[ii];
 628       distance.setValue(i, ii, newdist[ii]);
 629       distance.setValue(ii, i, newdist[ii]);
 630     }
 631   }
 632
 633   /**
 634    * DOCUMENT ME!
 635    *
 636    * @param i
 637    *          DOCUMENT ME!
 638    * @param j
 639    *          DOCUMENT ME!
 640    */
 641   public void findClusterNJDistance(int i, int j)
 642   {
 643
 644     // New distances from cluster to others
 645     double[] newdist = new double[noseqs];
 646
 647     for (int l = 0; l < noseqs; l++)
 648     {
 649       if ((l != i) && (l != j))
 650       {
 651         // newdist[l] = ((distance[i][l] + distance[j][l]) - distance[i][j]) /
 652         // 2;
 653         newdist[l] = (distance.getValue(i, l) + distance.getValue(j, l) - distance
 654                 .getValue(i, j)) / 2;
 655       }
 656       else
 657       {
 658         newdist[l] = 0;
 659       }
 660     }
 661
 662     for (int ii = 0; ii < noseqs; ii++)
 663     {
 664       // distance[i][ii] = newdist[ii];
 665       // distance[ii][i] = newdist[ii];
 666       distance.setValue(i, ii, newdist[ii]);
 667       distance.setValue(ii, i, newdist[ii]);
 668     }
 669   }
 670
 671   /**
 672    * DOCUMENT ME!
 673    *
 674    * @param i
 675    *          DOCUMENT ME!
 676    * @param j
 677    *          DOCUMENT ME!
 678    *
 679    * @return DOCUMENT ME!
 680    */
 681   public double findr(int i, int j)
 682   {
 683     double tmp = 1;
 684
 685     for (int k = 0; k < noseqs; k++)
 686     {
 687       if ((k != i) && (k != j) && (done[k] != 1))
 688       {
 689         // tmp = tmp + distance[i][k];
 690         tmp = tmp + distance.getValue(i, k);
 691       }
 692     }
 693
 694     if (noClus > 2)
 695     {
 696       tmp = tmp / (noClus - 2);
 697     }
 698
 699     return tmp;
 700   }
 701
 702   /**
 703    * DOCUMENT ME!
 704    *
 705    * @return DOCUMENT ME!
 706    */
 707   public double findMinNJDistance()
 708   {
 709     double min = Double.MAX_VALUE;
 710
 711     for (int i = 0; i < (noseqs - 1); i++)
 712     {
 713       for (int j = i + 1; j < noseqs; j++)
 714       {
 715         if ((done[i] != 1) && (done[j] != 1))
 716         {
 717           // float tmp = distance[i][j] - (findr(i, j) + findr(j, i));
 718           double tmp = distance.getValue(i, j)
 719                   - (findr(i, j) + findr(j, i));
 720
 721           if (tmp < min)
 722           {
 723             mini = i;
 724             minj = j;
 725
 726             min = tmp;
 727           }
 728         }
 729       }
 730     }
 731
 732     return min;
 733   }
 734
 735   /**
 736    * DOCUMENT ME!
 737    *
 738    * @return DOCUMENT ME!
 739    */
 740   public double findMinDistance()
 741   {
 742     double min = Double.MAX_VALUE;
 743
 744     for (int i = 0; i < (noseqs - 1); i++)
 745     {
 746       for (int j = i + 1; j < noseqs; j++)
 747       {
 748         if ((done[i] != 1) && (done[j] != 1))
 749         {
 750           // if (distance[i][j] < min)
 751           if (distance.getValue(i, j) < min)
 752           {
 753             mini = i;
 754             minj = j;
 755
 756             // min = distance[i][j];
 757             min = distance.getValue(i, j);
 758           }
 759         }
 760       }
 761     }
 762
 763     return min;
 764   }
 765
 766   /**
 767    * Calculate a distance matrix given the sequence input data and score model
 768    *
 769    * @return
 770    */
 771   public MatrixI findDistances(ScoreModelI scoreModel)
 772   {
 773     MatrixI result = null;
 774
 775     if (scoreModel == null)
 776     {
 777       // Resolve substitution model
 778       scoreModel = ScoreModels.getInstance().forName(pwtype);
 779       if (scoreModel == null)
 780       {
 781         scoreModel = ScoreModels.getInstance().forName("BLOSUM62");
 782       }
 783     }
 784     if (scoreModel instanceof DistanceScoreModelI)
 785     {
 786       result = ((DistanceScoreModelI) scoreModel).findDistances(seqData);
 787     }
 788     else if (scoreModel instanceof SimilarityScoreModelI)
 789     {
 790       /*
 791        * compute similarity and invert it to give a distance measure
 792        */
 793       result = ((SimilarityScoreModelI) scoreModel)
 794               .findSimilarities(seqData);
 795       result.reverseRange(true);
 796     }
 797     else
 798     {
 799       System.err
 800               .println("Unexpected type of score model, can't compute distances");
 801     }
 802
 803     return result;
 804   }
 805
 806   /**
 807    * DOCUMENT ME!
 808    */
 809   public void makeLeaves()
 810   {
 811     cluster = new Vector<Cluster>();
 812
 813     for (int i = 0; i < noseqs; i++)
 814     {
 815       SequenceNode sn = new SequenceNode();
 816
 817       sn.setElement(sequence[i]);
 818       sn.setName(sequence[i].getName());
 819       node.addElement(sn);
 820
 821       int[] value = new int[1];
 822       value[0] = i;
 823
 824       Cluster c = new Cluster(value);
 825       cluster.addElement(c);
 826     }
 827   }
 828
 829   /**
 830    * Search for leaf nodes below (or at) the given node
 831    *
 832    * @param nd
 833    *          root node to search from
 834    *
 835    * @return
 836    */
 837   public Vector<SequenceNode> findLeaves(SequenceNode nd)
 838   {
 839     Vector<SequenceNode> leaves = new Vector<SequenceNode>();
 840     findLeaves(nd, leaves);
 841     return leaves;
 842   }
 843
 844   /**
 845    * Search for leaf nodes.
 846    *
 847    * @param nd
 848    *          root node to search from
 849    * @param leaves
 850    *          Vector of leaves to add leaf node objects too.
 851    *
 852    * @return Vector of leaf nodes on binary tree
 853    */
 854   Vector<SequenceNode> findLeaves(SequenceNode nd,
 855           Vector<SequenceNode> leaves)
 856   {
 857     if (nd == null)
 858     {
 859       return leaves;
 860     }
 861
 862     if ((nd.left() == null) && (nd.right() == null)) // Interior node
 863     // detection
 864     {
 865       leaves.addElement(nd);
 866
 867       return leaves;
 868     }
 869     else
 870     {
 871       /*
 872        * TODO: Identify internal nodes... if (node.isSequenceLabel()) {
 873        * leaves.addElement(node); }
 874        */
 875       findLeaves((SequenceNode) nd.left(), leaves);
 876       findLeaves((SequenceNode) nd.right(), leaves);
 877     }
 878
 879     return leaves;
 880   }
 881
 882   /**
 883    * Find the leaf node with a particular ycount
 884    *
 885    * @param nd
 886    *          initial point on tree to search from
 887    * @param count
 888    *          value to search for
 889    *
 890    * @return null or the node with ycound=count
 891    */
 892   public Object findLeaf(SequenceNode nd, int count)
 893   {
 894     found = _findLeaf(nd, count);
 895
 896     return found;
 897   }
 898
 899   /*
 900    * #see findLeaf(SequenceNode node, count)
 901    */
 902   public Object _findLeaf(SequenceNode nd, int count)
 903   {
 904     if (nd == null)
 905     {
 906       return null;
 907     }
 908
 909     if (nd.ycount == count)
 910     {
 911       found = nd.element();
 912
 913       return found;
 914     }
 915     else
 916     {
 917       _findLeaf((SequenceNode) nd.left(), count);
 918       _findLeaf((SequenceNode) nd.right(), count);
 919     }
 920
 921     return found;
 922   }
 923
 924   /**
 925    * printNode is mainly for debugging purposes.
 926    *
 927    * @param nd
 928    *          SequenceNode
 929    */
 930   public void printNode(SequenceNode nd)
 931   {
 932     if (nd == null)
 933     {
 934       return;
 935     }
 936
 937     if ((nd.left() == null) && (nd.right() == null))
 938     {
 939       System.out.println("Leaf = " + ((SequenceI) nd.element()).getName());
 940       System.out.println("Dist " + nd.dist);
 941       System.out.println("Boot " + nd.getBootstrap());
 942     }
 943     else
 944     {
 945       System.out.println("Dist " + nd.dist);
 946       printNode((SequenceNode) nd.left());
 947       printNode((SequenceNode) nd.right());
 948     }
 949   }
 950
 951   /**
 952    * DOCUMENT ME!
 953    *
 954    * @param nd
 955    *          DOCUMENT ME!
 956    */
 957   public void findMaxDist(SequenceNode nd)
 958   {
 959     if (nd == null)
 960     {
 961       return;
 962     }
 963
 964     if ((nd.left() == null) && (nd.right() == null))
 965     {
 966       double dist = nd.dist;
 967
 968       if (dist > maxDistValue)
 969       {
 970         maxdist = nd;
 971         maxDistValue = dist;
 972       }
 973     }
 974     else
 975     {
 976       findMaxDist((SequenceNode) nd.left());
 977       findMaxDist((SequenceNode) nd.right());
 978     }
 979   }
 980
 981   /**
 982    * DOCUMENT ME!
 983    *
 984    * @return DOCUMENT ME!
 985    */
 986   public Vector<SequenceNode> getGroups()
 987   {
 988     return groups;
 989   }
 990
 991   /**
 992    * DOCUMENT ME!
 993    *
 994    * @return DOCUMENT ME!
 995    */
 996   public double getMaxHeight()
 997   {
 998     return maxheight;
 999   }
1000
1001   /**
1002    * DOCUMENT ME!
1003    *
1004    * @param nd
1005    *          DOCUMENT ME!
1006    * @param threshold
1007    *          DOCUMENT ME!
1008    */
1009   public void groupNodes(SequenceNode nd, float threshold)
1010   {
1011     if (nd == null)
1012     {
1013       return;
1014     }
1015
1016     if ((nd.height / maxheight) > threshold)
1017     {
1018       groups.addElement(nd);
1019     }
1020     else
1021     {
1022       groupNodes((SequenceNode) nd.left(), threshold);
1023       groupNodes((SequenceNode) nd.right(), threshold);
1024     }
1025   }
1026
1027   /**
1028    * DOCUMENT ME!
1029    *
1030    * @param nd
1031    *          DOCUMENT ME!
1032    *
1033    * @return DOCUMENT ME!
1034    */
1035   public double findHeight(SequenceNode nd)
1036   {
1037     if (nd == null)
1038     {
1039       return maxheight;
1040     }
1041
1042     if ((nd.left() == null) && (nd.right() == null))
1043     {
1044       nd.height = ((SequenceNode) nd.parent()).height + nd.dist;
1045
1046       if (nd.height > maxheight)
1047       {
1048         return nd.height;
1049       }
1050       else
1051       {
1052         return maxheight;
1053       }
1054     }
1055     else
1056     {
1057       if (nd.parent() != null)
1058       {
1059         nd.height = ((SequenceNode) nd.parent()).height + nd.dist;
1060       }
1061       else
1062       {
1063         maxheight = 0;
1064         nd.height = (float) 0.0;
1065       }
1066
1067       maxheight = findHeight((SequenceNode) (nd.left()));
1068       maxheight = findHeight((SequenceNode) (nd.right()));
1069     }
1070
1071     return maxheight;
1072   }
1073
1074   /**
1075    * DOCUMENT ME!
1076    *
1077    * @return DOCUMENT ME!
1078    */
1079   public SequenceNode reRoot()
1080   {
1081     if (maxdist != null)
1082     {
1083       ycount = 0;
1084
1085       double tmpdist = maxdist.dist;
1086
1087       // New top
1088       SequenceNode sn = new SequenceNode();
1089       sn.setParent(null);
1090
1091       // New right hand of top
1092       SequenceNode snr = (SequenceNode) maxdist.parent();
1093       changeDirection(snr, maxdist);
1094       System.out.println("Printing reversed tree");
1095       printN(snr);
1096       snr.dist = tmpdist / 2;
1097       maxdist.dist = tmpdist / 2;
1098
1099       snr.setParent(sn);
1100       maxdist.setParent(sn);
1101
1102       sn.setRight(snr);
1103       sn.setLeft(maxdist);
1104
1105       top = sn;
1106
1107       ycount = 0;
1108       reCount(top);
1109       findHeight(top);
1110     }
1111
1112     return top;
1113   }
1114
1115   /**
1116    *
1117    * @return true if original sequence data can be recovered
1118    */
1119   public boolean hasOriginalSequenceData()
1120   {
1121     return seqData != null;
1122   }
1123
1124   /**
1125    * Returns original alignment data used for calculation - or null where not
1126    * available.
1127    *
1128    * @return null or cut'n'pasteable alignment
1129    */
1130   public String printOriginalSequenceData(char gapChar)
1131   {
1132     if (seqData == null)
1133     {
1134       return null;
1135     }
1136
1137     StringBuffer sb = new StringBuffer();
1138     String[] seqdatas = seqData.getSequenceStrings(gapChar);
1139     for (int i = 0; i < seqdatas.length; i++)
1140     {
1141       sb.append(new jalview.util.Format("%-" + 15 + "s").form(sequence[i]
1142               .getName()));
1143       sb.append(" " + seqdatas[i] + "\n");
1144     }
1145     return sb.toString();
1146   }
1147
1148   /**
1149    * DOCUMENT ME!
1150    *
1151    * @param nd
1152    *          DOCUMENT ME!
1153    */
1154   public void printN(SequenceNode nd)
1155   {
1156     if (nd == null)
1157     {
1158       return;
1159     }
1160
1161     if ((nd.left() != null) && (nd.right() != null))
1162     {
1163       printN((SequenceNode) nd.left());
1164       printN((SequenceNode) nd.right());
1165     }
1166     else
1167     {
1168       System.out.println(" name = " + ((SequenceI) nd.element()).getName());
1169     }
1170
1171     System.out.println(" dist = " + nd.dist + " " + nd.count + " "
1172             + nd.height);
1173   }
1174
1175   /**
1176    * DOCUMENT ME!
1177    *
1178    * @param nd
1179    *          DOCUMENT ME!
1180    */
1181   public void reCount(SequenceNode nd)
1182   {
1183     ycount = 0;
1184     _lycount = 0;
1185     // _lylimit = this.node.size();
1186     _reCount(nd);
1187   }
1188
1189   private long _lycount = 0, _lylimit = 0;
1190
1191   /**
1192    * DOCUMENT ME!
1193    *
1194    * @param nd
1195    *          DOCUMENT ME!
1196    */
1197   public void _reCount(SequenceNode nd)
1198   {
1199     // if (_lycount<_lylimit)
1200     // {
1201     // System.err.println("Warning: depth of _recount greater than number of nodes.");
1202     // }
1203     if (nd == null)
1204     {
1205       return;
1206     }
1207     _lycount++;
1208
1209     if ((nd.left() != null) && (nd.right() != null))
1210     {
1211
1212       _reCount((SequenceNode) nd.left());
1213       _reCount((SequenceNode) nd.right());
1214
1215       SequenceNode l = (SequenceNode) nd.left();
1216       SequenceNode r = (SequenceNode) nd.right();
1217
1218       nd.count = l.count + r.count;
1219       nd.ycount = (l.ycount + r.ycount) / 2;
1220     }
1221     else
1222     {
1223       nd.count = 1;
1224       nd.ycount = ycount++;
1225     }
1226     _lycount--;
1227   }
1228
1229   /**
1230    * DOCUMENT ME!
1231    *
1232    * @param nd
1233    *          DOCUMENT ME!
1234    */
1235   public void swapNodes(SequenceNode nd)
1236   {
1237     if (nd == null)
1238     {
1239       return;
1240     }
1241
1242     SequenceNode tmp = (SequenceNode) nd.left();
1243
1244     nd.setLeft(nd.right());
1245     nd.setRight(tmp);
1246   }
1247
1248   /**
1249    * DOCUMENT ME!
1250    *
1251    * @param nd
1252    *          DOCUMENT ME!
1253    * @param dir
1254    *          DOCUMENT ME!
1255    */
1256   public void changeDirection(SequenceNode nd, SequenceNode dir)
1257   {
1258     if (nd == null)
1259     {
1260       return;
1261     }
1262
1263     if (nd.parent() != top)
1264     {
1265       changeDirection((SequenceNode) nd.parent(), nd);
1266
1267       SequenceNode tmp = (SequenceNode) nd.parent();
1268
1269       if (dir == nd.left())
1270       {
1271         nd.setParent(dir);
1272         nd.setLeft(tmp);
1273       }
1274       else if (dir == nd.right())
1275       {
1276         nd.setParent(dir);
1277         nd.setRight(tmp);
1278       }
1279     }
1280     else
1281     {
1282       if (dir == nd.left())
1283       {
1284         nd.setParent(nd.left());
1285
1286         if (top.left() == nd)
1287         {
1288           nd.setRight(top.right());
1289         }
1290         else
1291         {
1292           nd.setRight(top.left());
1293         }
1294       }
1295       else
1296       {
1297         nd.setParent(nd.right());
1298
1299         if (top.left() == nd)
1300         {
1301           nd.setLeft(top.right());
1302         }
1303         else
1304         {
1305           nd.setLeft(top.left());
1306         }
1307       }
1308     }
1309   }
1310
1311   /**
1312    * DOCUMENT ME!
1313    *
1314    * @return DOCUMENT ME!
1315    */
1316   public SequenceNode getMaxDist()
1317   {
1318     return maxdist;
1319   }
1320
1321   /**
1322    * DOCUMENT ME!
1323    *
1324    * @return DOCUMENT ME!
1325    */
1326   public SequenceNode getTopNode()
1327   {
1328     return top;
1329   }
1330
1331   /**
1332    *
1333    * @return true if tree has real distances
1334    */
1335   public boolean isHasDistances()
1336   {
1337     return hasDistances;
1338   }
1339
1340   /**
1341    *
1342    * @return true if tree has real bootstrap values
1343    */
1344   public boolean isHasBootstrap()
1345   {
1346     return hasBootstrap;
1347   }
1348
1349   public boolean isHasRootDistance()
1350   {
1351     return hasRootDistance;
1352   }
1353
1354   /**
1355    * apply the given transform to all the nodes in the tree.
1356    *
1357    * @param nodeTransformI
1358    */
1359   public void applyToNodes(NodeTransformI nodeTransformI)
1360   {
1361     for (Enumeration<SequenceNode> nodes = node.elements(); nodes
1362             .hasMoreElements(); nodeTransformI.transform(nodes
1363             .nextElement()))
1364     {
1365       ;
1366     }
1367   }
1368 }
1369
1370 /**
1371  * DOCUMENT ME!
1372  *
1373  * @author $author$
1374  * @version $Revision$
1375  */
1376 // TODO what does this class have that int[] doesn't have already?
1377 class Cluster
1378 {
1379   int[] value;
1380
1381   /**
1382    * Creates a new Cluster object.
1383    *
1384    * @param value
1385    *          DOCUMENT ME!
1386    */
1387   public Cluster(int[] value)
1388   {
1389     this.value = value;
1390   }
1391 }