src/jalview/analysis/Conservation.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 java.util.Locale;
  24
  25 import jalview.analysis.scoremodels.ScoreMatrix;
  26 import jalview.analysis.scoremodels.ScoreModels;
  27 import jalview.api.analysis.ScoreModelI;
  28 import jalview.datamodel.AlignmentAnnotation;
  29 import jalview.datamodel.Annotation;
  30 import jalview.datamodel.ResidueCount;
  31 import jalview.datamodel.ResidueCount.SymbolCounts;
  32 import jalview.datamodel.Sequence;
  33 import jalview.datamodel.SequenceI;
  34 import jalview.schemes.ResidueProperties;
  35 import jalview.util.Comparison;
  36 import jalview.util.Format;
  37
  38 import java.awt.Color;
  39 import java.util.List;
  40 import java.util.Map;
  41 import java.util.Map.Entry;
  42 import java.util.SortedMap;
  43 import java.util.TreeMap;
  44 import java.util.Vector;
  45
  46 /**
  47  * Calculates conservation values for a given set of sequences
  48  */
  49 public class Conservation
  50 {
  51   /*
  52    * need to have a minimum of 3% of sequences with a residue
  53    * for it to be included in the conservation calculation
  54    */
  55   private static final int THRESHOLD_PERCENT = 3;
  56
  57   private static final int TOUPPERCASE = 'a' - 'A';
  58
  59   private static final int GAP_INDEX = -1;
  60
  61   private static final Format FORMAT_3DP = new Format("%2.5f");
  62
  63   SequenceI[] sequences;
  64
  65   int start;
  66
  67   int end;
  68
  69   /*
  70    * a list whose i'th element is an array whose first entry is the checksum
  71    * of the i'th sequence, followed by residues encoded to score matrix index
  72    */
  73   Vector<int[]> seqNums;
  74
  75   int maxLength = 0; // used by quality calcs
  76
  77   boolean seqNumsChanged = false; // updated after any change via calcSeqNum;
  78
  79   /*
  80    * a map per column with {property, conservation} where conservation value is
  81    * 1 (property is conserved), 0 (absence of property is conserved) or -1
  82    * (property is not conserved i.e. column has residues with and without it)
  83    */
  84   Map<String, Integer>[] total;
  85
  86   /*
  87    * if true then conservation calculation will map all symbols to canonical aa
  88    * numbering rather than consider conservation of that symbol
  89    */
  90   boolean canonicaliseAa = true;
  91
  92   private Vector<Double> quality;
  93
  94   private double qualityMinimum;
  95
  96   private double qualityMaximum;
  97
  98   private Sequence consSequence;
  99
 100   /*
 101    * percentage of residues in a column to qualify for counting conservation
 102    */
 103   private int threshold;
 104
 105   private String name = "";
 106
 107   /*
 108    * an array, for each column, of counts of symbols (by score matrix index)
 109    */
 110   private int[][] cons2;
 111
 112   /*
 113    * gap counts for each column
 114    */
 115   private int[] cons2GapCounts;
 116
 117   private String[] consSymbs;
 118
 119   /**
 120    * Constructor using default threshold of 3%
 121    *
 122    * @param name
 123    *          Name of conservation
 124    * @param sequences
 125    *          sequences to be used in calculation
 126    * @param start
 127    *          start residue position
 128    * @param end
 129    *          end residue position
 130    */
 131   public Conservation(String name, List<SequenceI> sequences, int start,
 132           int end)
 133   {
 134     this(name, THRESHOLD_PERCENT, sequences, start, end);
 135   }
 136
 137   /**
 138    * Constructor
 139    *
 140    * @param name
 141    *          Name of conservation
 142    * @param threshold
 143    *          percentage of sequences at or below which property conservation is
 144    *          ignored
 145    * @param sequences
 146    *          sequences to be used in calculation
 147    * @param start
 148    *          start column position
 149    * @param end
 150    *          end column position
 151    */
 152   public Conservation(String name, int threshold, List<SequenceI> sequences,
 153           int start, int end)
 154   {
 155     this.name = name;
 156     this.threshold = threshold;
 157     this.start = start;
 158     this.end = end;
 159
 160     maxLength = end - start + 1; // default width includes bounds of
 161     // calculation
 162
 163     int s, sSize = sequences.size();
 164     SequenceI[] sarray = new SequenceI[sSize];
 165     this.sequences = sarray;
 166     try
 167     {
 168       for (s = 0; s < sSize; s++)
 169       {
 170         sarray[s] = sequences.get(s);
 171         if (sarray[s].getLength() > maxLength)
 172         {
 173           maxLength = sarray[s].getLength();
 174         }
 175       }
 176     } catch (ArrayIndexOutOfBoundsException ex)
 177     {
 178       // bail - another thread has modified the sequence array, so the current
 179       // calculation is probably invalid.
 180       this.sequences = new SequenceI[0];
 181       maxLength = 0;
 182     }
 183   }
 184
 185   /**
 186    * Translate sequence i into score matrix indices and store it in the i'th
 187    * position of the seqNums array.
 188    *
 189    * @param i
 190    * @param sm
 191    */
 192   private void calcSeqNum(int i, ScoreMatrix sm)
 193   {
 194     int sSize = sequences.length;
 195
 196     if ((i > -1) && (i < sSize))
 197     {
 198       String sq = sequences[i].getSequenceAsString();
 199
 200       if (seqNums.size() <= i)
 201       {
 202         seqNums.addElement(new int[sq.length() + 1]);
 203       }
 204
 205       /*
 206        * the first entry in the array is the sequence's hashcode,
 207        * following entries are matrix indices of sequence characters
 208        */
 209       if (sq.hashCode() != seqNums.elementAt(i)[0])
 210       {
 211         int j;
 212         int len;
 213         seqNumsChanged = true;
 214         len = sq.length();
 215
 216         if (maxLength < len)
 217         {
 218           maxLength = len;
 219         }
 220
 221         int[] sqnum = new int[len + 1]; // better to always make a new array -
 222         // sequence can change its length
 223         sqnum[0] = sq.hashCode();
 224
 225         for (j = 1; j <= len; j++)
 226         {
 227           // sqnum[j] = ResidueProperties.aaIndex[sq.charAt(j - 1)];
 228           char residue = sq.charAt(j - 1);
 229           if (Comparison.isGap(residue))
 230           {
 231             sqnum[j] = GAP_INDEX;
 232           }
 233           else
 234           {
 235             sqnum[j] = sm.getMatrixIndex(residue);
 236             if (sqnum[j] == -1)
 237             {
 238               sqnum[j] = GAP_INDEX;
 239             }
 240           }
 241         }
 242
 243         seqNums.setElementAt(sqnum, i);
 244       }
 245       else
 246       {
 247         jalview.bin.Console.outPrintln("SEQUENCE HAS BEEN DELETED!!!");
 248       }
 249     }
 250     else
 251     {
 252       // JBPNote INFO level debug
 253       jalview.bin.Console.errPrintln(
 254               "ERROR: calcSeqNum called with out of range sequence index for Alignment\n");
 255     }
 256   }
 257
 258   /**
 259    * Calculates the conservation values for given set of sequences
 260    */
 261   public void calculate()
 262   {
 263     int height = sequences.length;
 264
 265     total = new Map[maxLength];
 266
 267     for (int column = start; column <= end; column++)
 268     {
 269       ResidueCount values = countResidues(column);
 270
 271       /*
 272        * percentage count at or below which we ignore residues
 273        */
 274       int thresh = (threshold * height) / 100;
 275
 276       /*
 277        * check observed residues in column and record whether each
 278        * physico-chemical property is conserved (+1), absence conserved (0),
 279        * or not conserved (-1)
 280        * Using TreeMap means properties are displayed in alphabetical order
 281        */
 282       SortedMap<String, Integer> resultHash = new TreeMap<>();
 283       SymbolCounts symbolCounts = values.getSymbolCounts();
 284       char[] symbols = symbolCounts.symbols;
 285       int[] counts = symbolCounts.values;
 286       for (int j = 0; j < symbols.length; j++)
 287       {
 288         char c = symbols[j];
 289         if (counts[j] > thresh)
 290         {
 291           recordConservation(resultHash, String.valueOf(c));
 292         }
 293       }
 294       if (values.getGapCount() > thresh)
 295       {
 296         recordConservation(resultHash, "-");
 297       }
 298
 299       if (total.length > 0)
 300       {
 301         total[column - start] = resultHash;
 302       }
 303     }
 304   }
 305
 306   /**
 307    * Updates the conservation results for an observed residue
 308    *
 309    * @param resultMap
 310    *          a map of {property, conservation} where conservation value is +1
 311    *          (all residues have the property), 0 (no residue has the property)
 312    *          or -1 (some do, some don't)
 313    * @param res
 314    */
 315   protected static void recordConservation(Map<String, Integer> resultMap,
 316           String res)
 317   {
 318     res = res.toUpperCase(Locale.ROOT);
 319     for (Entry<String, Map<String, Integer>> property : ResidueProperties.propHash
 320             .entrySet())
 321     {
 322       String propertyName = property.getKey();
 323       Integer residuePropertyValue = property.getValue().get(res);
 324
 325       if (!resultMap.containsKey(propertyName))
 326       {
 327         /*
 328          * first time we've seen this residue - note whether it has this property
 329          */
 330         if (residuePropertyValue != null)
 331         {
 332           resultMap.put(propertyName, residuePropertyValue);
 333         }
 334         else
 335         {
 336           /*
 337            * unrecognised residue - use default value for property
 338            */
 339           resultMap.put(propertyName, property.getValue().get("-"));
 340         }
 341       }
 342       else
 343       {
 344         Integer currentResult = resultMap.get(propertyName);
 345         if (currentResult.intValue() != -1
 346                 && !currentResult.equals(residuePropertyValue))
 347         {
 348           /*
 349            * property is unconserved - residues seen both with and without it
 350            */
 351           resultMap.put(propertyName, Integer.valueOf(-1));
 352         }
 353       }
 354     }
 355   }
 356
 357   /**
 358    * Counts residues (upper-cased) and gaps in the given column
 359    *
 360    * @param column
 361    * @return
 362    */
 363   protected ResidueCount countResidues(int column)
 364   {
 365     ResidueCount values = new ResidueCount(false);
 366
 367     for (int row = 0; row < sequences.length; row++)
 368     {
 369       if (sequences[row].getLength() > column)
 370       {
 371         char c = sequences[row].getCharAt(column);
 372         if (canonicaliseAa)
 373         {
 374           int index = ResidueProperties.aaIndex[c];
 375           c = index > 20 ? '-' : ResidueProperties.aa[index].charAt(0);
 376         }
 377         else
 378         {
 379           c = toUpperCase(c);
 380         }
 381         if (Comparison.isGap(c))
 382         {
 383           values.addGap();
 384         }
 385         else
 386         {
 387           values.add(c);
 388         }
 389       }
 390       else
 391       {
 392         values.addGap();
 393       }
 394     }
 395     return values;
 396   }
 397
 398   /**
 399    * Counts conservation and gaps for a column of the alignment
 400    *
 401    * @return { 1 if fully conserved, else 0, gap count }
 402    */
 403   public int[] countConservationAndGaps(int column)
 404   {
 405     int gapCount = 0;
 406     boolean fullyConserved = true;
 407     int iSize = sequences.length;
 408
 409     if (iSize == 0)
 410     {
 411       return new int[] { 0, 0 };
 412     }
 413
 414     char lastRes = '0';
 415     for (int i = 0; i < iSize; i++)
 416     {
 417       if (column >= sequences[i].getLength())
 418       {
 419         gapCount++;
 420         continue;
 421       }
 422
 423       char c = sequences[i].getCharAt(column); // gaps do not have upper/lower
 424                                                // case
 425
 426       if (Comparison.isGap((c)))
 427       {
 428         gapCount++;
 429       }
 430       else
 431       {
 432         c = toUpperCase(c);
 433         if (lastRes == '0')
 434         {
 435           lastRes = c;
 436         }
 437         if (c != lastRes)
 438         {
 439           fullyConserved = false;
 440         }
 441       }
 442     }
 443
 444     int[] r = new int[] { fullyConserved ? 1 : 0, gapCount };
 445     return r;
 446   }
 447
 448   /**
 449    * Returns the upper-cased character if between 'a' and 'z', else the
 450    * unchanged value
 451    *
 452    * @param c
 453    * @return
 454    */
 455   char toUpperCase(char c)
 456   {
 457     if ('a' <= c && c <= 'z')
 458     {
 459       c -= TOUPPERCASE;
 460     }
 461     return c;
 462   }
 463
 464   /**
 465    * Calculates the conservation sequence
 466    *
 467    * @param positiveOnly
 468    *          if true, calculate positive conservation; else calculate both
 469    *          positive and negative conservation
 470    * @param maxPercentageGaps
 471    *          the percentage of gaps in a column, at or above which no
 472    *          conservation is asserted
 473    */
 474   public void verdict(boolean positiveOnly, float maxPercentageGaps)
 475   {
 476     // TODO call this at the end of calculate(), should not be a public method
 477
 478     StringBuilder consString = new StringBuilder(end);
 479
 480     // NOTE THIS SHOULD CHECK IF THE CONSEQUENCE ALREADY
 481     // EXISTS AND NOT OVERWRITE WITH '-', BUT THIS CASE
 482     // DOES NOT EXIST IN JALVIEW 2.1.2
 483     for (int i = 0; i < start; i++)
 484     {
 485       consString.append('-');
 486     }
 487     consSymbs = new String[end - start + 1];
 488     for (int i = start; i <= end; i++)
 489     {
 490       int[] gapcons = countConservationAndGaps(i);
 491       boolean fullyConserved = gapcons[0] == 1;
 492       int totGaps = gapcons[1];
 493       float pgaps = (totGaps * 100f) / sequences.length;
 494
 495       if (maxPercentageGaps > pgaps)
 496       {
 497         Map<String, Integer> resultHash = total[i - start];
 498         int count = 0;
 499         StringBuilder positives = new StringBuilder(64);
 500         StringBuilder negatives = new StringBuilder(32);
 501         for (String type : resultHash.keySet())
 502         {
 503           int result = resultHash.get(type).intValue();
 504           if (result == -1)
 505           {
 506             /*
 507              * not conserved (present or absent)
 508              */
 509             continue;
 510           }
 511           count++;
 512           if (result == 1)
 513           {
 514             /*
 515              * positively conserved property (all residues have it)
 516              */
 517             positives.append(positives.length() == 0 ? "" : " ");
 518             positives.append(type);
 519           }
 520           if (result == 0 && !positiveOnly)
 521           {
 522             /*
 523              * absense of property is conserved (all residues lack it)
 524              */
 525             negatives.append(negatives.length() == 0 ? "" : " ");
 526             negatives.append("!").append(type);
 527           }
 528         }
 529         if (negatives.length() > 0)
 530         {
 531           positives.append(" ").append(negatives);
 532         }
 533         consSymbs[i - start] = positives.toString();
 534
 535         if (count < 10)
 536         {
 537           consString.append(count); // Conserved props!=Identity
 538         }
 539         else
 540         {
 541           consString.append(fullyConserved ? "*" : "+");
 542         }
 543       }
 544       else
 545       {
 546         consString.append('-');
 547       }
 548     }
 549
 550     consSequence = new Sequence(name, consString.toString(), start, end);
 551   }
 552
 553   /**
 554    *
 555    *
 556    * @return Conservation sequence
 557    */
 558   public SequenceI getConsSequence()
 559   {
 560     return consSequence;
 561   }
 562
 563   // From Alignment.java in jalview118
 564   public void findQuality()
 565   {
 566     findQuality(0, maxLength - 1, getQualitySubstMat());
 567   }
 568
 569   public ScoreMatrix getQualitySubstMat()
 570   {
 571     return sm;
 572   }
 573
 574   /**
 575    * default for quality calc is BLOSUM
 576    */
 577   private ScoreMatrix sm = ScoreModels.getInstance().getBlosum62();
 578
 579   /**
 580    * set the matrix used for calculating quality scores
 581    *
 582    * @param scoreMatrix
 583    *          - must be a valid score matrix able to evaluate all amino acid
 584    *          substitutions
 585    */
 586   public void setQualitySubstMat(ScoreMatrix scoreMatrix)
 587   {
 588     if (scoreMatrix != null && scoreMatrix.isProtein())
 589     {
 590       sm = scoreMatrix;
 591     }
 592   }
 593
 594   /**
 595    * DOCUMENT ME!
 596    *
 597    * @param sm
 598    */
 599   private void percentIdentity(ScoreMatrix sm)
 600   {
 601     seqNums = new Vector<>();
 602     int i = 0, iSize = sequences.length;
 603     // Do we need to calculate this again?
 604     for (i = 0; i < iSize; i++)
 605     {
 606       calcSeqNum(i, sm);
 607     }
 608
 609     if ((cons2 == null) || seqNumsChanged)
 610     {
 611       // FIXME remove magic number 24 without changing calc
 612       // sm.getSize() returns 25 so doesn't quite do it...
 613       cons2 = new int[maxLength][24];
 614       cons2GapCounts = new int[maxLength];
 615
 616       int j = 0;
 617
 618       while (j < sequences.length)
 619       {
 620         int[] sqnum = seqNums.elementAt(j);
 621
 622         for (i = 1; i < sqnum.length; i++)
 623         {
 624           int index = sqnum[i];
 625           if (index == GAP_INDEX)
 626           {
 627             cons2GapCounts[i - 1]++;
 628           }
 629           else
 630           {
 631             cons2[i - 1][index]++;
 632           }
 633         }
 634
 635         // TODO should this start from sqnum.length?
 636         for (i = sqnum.length - 1; i < maxLength; i++)
 637         {
 638           cons2GapCounts[i]++;
 639         }
 640         j++;
 641       }
 642     }
 643   }
 644
 645   /**
 646    * Calculates the quality of the set of sequences over the given inclusive
 647    * column range, using the specified substitution score matrix
 648    *
 649    * @param startCol
 650    * @param endCol
 651    * @param scoreMatrix
 652    */
 653   protected void findQuality(int startCol, int endCol,
 654           ScoreMatrix scoreMatrix)
 655   {
 656     quality = new Vector<>();
 657
 658     double max = -Double.MAX_VALUE;
 659     float[][] scores = scoreMatrix.getMatrix();
 660
 661     percentIdentity(scoreMatrix);
 662
 663     int size = seqNums.size();
 664     int[] lengths = new int[size];
 665
 666     for (int l = 0; l < size; l++)
 667     {
 668       lengths[l] = seqNums.elementAt(l).length - 1;
 669     }
 670
 671     final int symbolCount = scoreMatrix.getSize();
 672
 673     for (int j = startCol; j <= endCol; j++)
 674     {
 675       double bigtot = 0;
 676
 677       // First Xr = depends on column only
 678       double[] x = new double[symbolCount];
 679
 680       for (int ii = 0; ii < symbolCount; ii++)
 681       {
 682         x[ii] = 0;
 683
 684         /*
 685          * todo JAL-728 currently assuming last symbol in matrix is * for gap
 686          * (which we ignore as counted separately); true for BLOSUM62 but may
 687          * not be once alternative matrices are supported
 688          */
 689         for (int i2 = 0; i2 < symbolCount - 1; i2++)
 690         {
 691           x[ii] += (((double) cons2[j][i2] * scores[ii][i2]) + 4D);
 692         }
 693         x[ii] += 4D + cons2GapCounts[j] * scoreMatrix.getMinimumScore();
 694
 695         x[ii] /= size;
 696       }
 697
 698       // Now calculate D for each position and sum
 699       for (int k = 0; k < size; k++)
 700       {
 701         double tot = 0;
 702         double[] xx = new double[symbolCount];
 703         // sequence character index, or implied gap if sequence too short
 704         int seqNum = (j < lengths[k]) ? seqNums.elementAt(k)[j + 1]
 705                 : GAP_INDEX;
 706
 707         for (int i = 0; i < symbolCount - 1; i++)
 708         {
 709           double sr = 4D;
 710           if (seqNum == GAP_INDEX)
 711           {
 712             sr += scoreMatrix.getMinimumScore();
 713           }
 714           else
 715           {
 716             sr += scores[i][seqNum];
 717           }
 718
 719           xx[i] = x[i] - sr;
 720
 721           tot += (xx[i] * xx[i]);
 722         }
 723
 724         bigtot += Math.sqrt(tot);
 725       }
 726
 727       max = Math.max(max, bigtot);
 728
 729       quality.addElement(Double.valueOf(bigtot));
 730     }
 731
 732     double newmax = -Double.MAX_VALUE;
 733
 734     for (int j = startCol; j <= endCol; j++)
 735     {
 736       double tmp = quality.elementAt(j).doubleValue();
 737       // tmp = ((max - tmp) * (size - cons2[j][23])) / size;
 738       tmp = ((max - tmp) * (size - cons2GapCounts[j])) / size;
 739
 740       // jalview.bin.Console.outPrintln(tmp+ " " + j);
 741       quality.setElementAt(Double.valueOf(tmp), j);
 742
 743       if (tmp > newmax)
 744       {
 745         newmax = tmp;
 746       }
 747     }
 748
 749     qualityMinimum = 0D;
 750     qualityMaximum = newmax;
 751   }
 752
 753   /**
 754    * Complete the given consensus and quuality annotation rows. Note: currently
 755    * this method will reallocate the given annotation row if it is different to
 756    * the calculated width, otherwise will leave its length unchanged.
 757    *
 758    * @param conservation
 759    *          conservation annotation row
 760    * @param quality2
 761    *          (optional - may be null)
 762    * @param istart
 763    *          first column for conservation
 764    * @param alWidth
 765    *          extent of conservation
 766    */
 767   public void completeAnnotations(AlignmentAnnotation conservation,
 768           AlignmentAnnotation quality2, int istart, int alWidth)
 769   {
 770     SequenceI cons = getConsSequence();
 771
 772     /*
 773      * colour scale for Conservation and Quality;
 774      */
 775     float minR = 0.3f;
 776     float minG = 0.0f;
 777     float minB = 0f;
 778     float maxR = 1.0f - minR;
 779     float maxG = 0.9f - minG;
 780     float maxB = 0f - minB;
 781
 782     float min = 0f;
 783     float max = 11f;
 784     float qmin = 0f;
 785     float qmax = 0f;
 786
 787     if (conservation != null && conservation.annotations != null
 788             && conservation.annotations.length != alWidth)
 789     {
 790       conservation.annotations = new Annotation[alWidth];
 791     }
 792
 793     if (quality2 != null)
 794     {
 795       quality2.graphMax = (float) qualityMaximum;
 796       if (quality2.annotations != null
 797               && quality2.annotations.length != alWidth)
 798       {
 799         quality2.annotations = new Annotation[alWidth];
 800       }
 801       qmin = (float) qualityMinimum;
 802       qmax = (float) qualityMaximum;
 803     }
 804
 805     for (int i = istart; i < alWidth; i++)
 806     {
 807       float value = 0;
 808
 809       char c = cons.getCharAt(i);
 810
 811       if (Character.isDigit(c))
 812       {
 813         value = c - '0';
 814       }
 815       else if (c == '*')
 816       {
 817         value = 11;
 818       }
 819       else if (c == '+')
 820       {
 821         value = 10;
 822       }
 823
 824       if (conservation != null)
 825       {
 826         float vprop = value - min;
 827         vprop /= max;
 828         int consp = i - start;
 829         String conssym = (value > 0 && consp > -1
 830                 && consp < consSymbs.length) ? consSymbs[consp] : "";
 831         conservation.annotations[i] = new Annotation(String.valueOf(c),
 832                 conssym, ' ', value, new Color(minR + (maxR * vprop),
 833                         minG + (maxG * vprop), minB + (maxB * vprop)));
 834       }
 835
 836       // Quality calc
 837       if (quality2 != null)
 838       {
 839         value = quality.elementAt(i).floatValue();
 840         float vprop = value - qmin;
 841         vprop /= qmax;
 842         String description = FORMAT_3DP.form(value);
 843         quality2.annotations[i] = new Annotation(" ", description, ' ',
 844                 value, new Color(minR + (maxR * vprop),
 845                         minG + (maxG * vprop), minB + (maxB * vprop)));
 846       }
 847     }
 848   }
 849
 850   /**
 851    * construct and call the calculation methods on a new Conservation object
 852    *
 853    * @param name
 854    *          - name of conservation
 855    * @param seqs
 856    * @param start
 857    *          first column in calculation window
 858    * @param end
 859    *          last column in calculation window
 860    * @param positiveOnly
 861    *          calculate positive (true) or positive and negative (false)
 862    *          conservation
 863    * @param maxPercentGaps
 864    *          percentage of gaps tolerated in column
 865    * @param calcQuality
 866    *          flag indicating if alignment quality should be calculated
 867    * @return Conservation object ready for use in visualization
 868    */
 869   public static Conservation calculateConservation(String name,
 870           List<SequenceI> seqs, int start, int end, boolean positiveOnly,
 871           int maxPercentGaps, boolean calcQuality)
 872   {
 873     return calculateConservation(name, seqs, start, end, positiveOnly,
 874             maxPercentGaps, calcQuality, null);
 875   }
 876
 877   /**
 878    * construct and call the calculation methods on a new Conservation object
 879    *
 880    * @param name
 881    *          - name of conservation
 882    * @param seqs
 883    * @param start
 884    *          first column in calculation window
 885    * @param end
 886    *          last column in calculation window
 887    * @param positiveOnly
 888    *          calculate positive (true) or positive and negative (false)
 889    *          conservation
 890    * @param maxPercentGaps
 891    *          percentage of gaps tolerated in column
 892    * @param calcQuality
 893    *          flag indicating if alignment quality should be calculated
 894    * @param scoreMatrix
 895    *          - null or peptide score matrix used for quality calculation
 896    * @return Conservation object ready for use in visualization
 897    */
 898   public static Conservation calculateConservation(String name,
 899           List<SequenceI> seqs, int start, int end, boolean positiveOnly,
 900           int maxPercentGaps, boolean calcQuality, ScoreMatrix scoreMatrix)
 901   {
 902     Conservation cons = new Conservation(name, seqs, start, end);
 903     if (scoreMatrix != null)
 904     {
 905       cons.setQualitySubstMat(scoreMatrix);
 906     }
 907     cons.calculate();
 908     cons.verdict(positiveOnly, maxPercentGaps);
 909
 910     if (calcQuality)
 911     {
 912       cons.findQuality();
 913     }
 914
 915     return cons;
 916   }
 917
 918   /**
 919    * Returns the computed tooltip (annotation description) for a given column.
 920    * The tip is empty if the conservation score is zero, otherwise holds the
 921    * conserved properties (and, optionally, properties whose absence is
 922    * conserved).
 923    *
 924    * @param column
 925    * @return
 926    */
 927   String getTooltip(int column)
 928   {
 929     SequenceI cons = getConsSequence();
 930     char val = column < cons.getLength() ? cons.getCharAt(column) : '-';
 931     boolean hasConservation = val != '-' && val != '0';
 932     int consp = column - start;
 933     String tip = (hasConservation && consp > -1 && consp < consSymbs.length)
 934             ? consSymbs[consp]
 935             : "";
 936     return tip;
 937   }
 938 }