*/
package jalview.analysis;
+import java.util.Locale;
+
import jalview.analysis.scoremodels.ScoreMatrix;
import jalview.analysis.scoremodels.ScoreModels;
import jalview.datamodel.AlignmentAnnotation;
import jalview.datamodel.SequenceI;
import jalview.schemes.ResidueProperties;
import jalview.util.Comparison;
+import jalview.util.Format;
import java.awt.Color;
import java.util.List;
private static final int TOUPPERCASE = 'a' - 'A';
+ private static final int GAP_INDEX = -1;
+
+ private static final Format FORMAT_3DP = new Format("%2.5f");
+
SequenceI[] sequences;
int start;
int end;
- Vector<int[]> seqNums; // vector of int vectors where first is sequence
- // checksum
+ /*
+ * a list whose i'th element is an array whose first entry is the checksum
+ * of the i'th sequence, followed by residues encoded to score matrix index
+ */
+ Vector<int[]> seqNums;
int maxLength = 0; // used by quality calcs
*/
Map<String, Integer>[] total;
- boolean canonicaliseAa = true; // if true then conservation calculation will
-
- // map all symbols to canonical aa numbering
- // rather than consider conservation of that
- // symbol
+ /*
+ * if true then conservation calculation will map all symbols to canonical aa
+ * numbering rather than consider conservation of that symbol
+ */
+ boolean canonicaliseAa = true;
- /** Stores calculated quality values */
private Vector<Double> quality;
- /** Stores maximum and minimum values of quality values */
- private double[] qualityRange = new double[2];
+ private double qualityMinimum;
+
+ private double qualityMaximum;
private Sequence consSequence;
private String name = "";
+ /*
+ * an array, for each column, of counts of symbols (by score matrix index)
+ */
private int[][] cons2;
+ /*
+ * gap counts for each column
+ */
+ private int[] cons2GapCounts;
+
private String[] consSymbs;
/**
* @param end
* end column position
*/
- public Conservation(String name, int threshold,
- List<SequenceI> sequences, int start, int end)
+ public Conservation(String name, int threshold, List<SequenceI> sequences,
+ int start, int end)
{
this.name = name;
this.threshold = threshold;
*/
private void calcSeqNum(int i, ScoreMatrix sm)
{
- int gapIndex = sm.getGapIndex();
int sSize = sequences.length;
if ((i > -1) && (i < sSize))
for (j = 1; j <= len; j++)
{
// sqnum[j] = ResidueProperties.aaIndex[sq.charAt(j - 1)];
- sqnum[j] = sm.getMatrixIndex(sq.charAt(j - 1));
- if (sqnum[j] == -1)
+ char residue = sq.charAt(j - 1);
+ if (Comparison.isGap(residue))
+ {
+ sqnum[j] = GAP_INDEX;
+ }
+ else
{
- sqnum[j] = gapIndex;
+ sqnum[j] = sm.getMatrixIndex(residue);
+ if (sqnum[j] == -1)
+ {
+ sqnum[j] = GAP_INDEX;
+ }
}
}
else
{
// JBPNote INFO level debug
- System.err
- .println("ERROR: calcSeqNum called with out of range sequence index for Alignment\n");
+ System.err.println(
+ "ERROR: calcSeqNum called with out of range sequence index for Alignment\n");
}
}
* or not conserved (-1)
* Using TreeMap means properties are displayed in alphabetical order
*/
- SortedMap<String, Integer> resultHash = new TreeMap<String, Integer>();
+ SortedMap<String, Integer> resultHash = new TreeMap<>();
SymbolCounts symbolCounts = values.getSymbolCounts();
char[] symbols = symbolCounts.symbols;
int[] counts = symbolCounts.values;
protected static void recordConservation(Map<String, Integer> resultMap,
String res)
{
- res = res.toUpperCase();
+ res = res.toUpperCase(Locale.ROOT);
for (Entry<String, Map<String, Integer>> property : ResidueProperties.propHash
.entrySet())
{
continue;
}
- char c = sequences[i].getCharAt(column); // gaps do not have upper/lower case
+ char c = sequences[i].getCharAt(column); // gaps do not have upper/lower
+ // case
if (Comparison.isGap((c)))
{
// From Alignment.java in jalview118
public void findQuality()
{
- findQuality(0, maxLength - 1);
+ findQuality(0, maxLength - 1, ScoreModels.getInstance().getBlosum62());
}
/**
*/
private void percentIdentity(ScoreMatrix sm)
{
- seqNums = new Vector<int[]>();
+ seqNums = new Vector<>();
int i = 0, iSize = sequences.length;
// Do we need to calculate this again?
for (i = 0; i < iSize; i++)
calcSeqNum(i, sm);
}
- int gapIndex = sm.getGapIndex();
-
if ((cons2 == null) || seqNumsChanged)
{
// FIXME remove magic number 24 without changing calc
// sm.getSize() returns 25 so doesn't quite do it...
cons2 = new int[maxLength][24];
+ cons2GapCounts = new int[maxLength];
int j = 0;
for (i = 1; i < sqnum.length; i++)
{
- cons2[i - 1][sqnum[i]]++;
+ int index = sqnum[i];
+ if (index == GAP_INDEX)
+ {
+ cons2GapCounts[i - 1]++;
+ }
+ else
+ {
+ cons2[i - 1][index]++;
+ }
}
+ // TODO should this start from sqnum.length?
for (i = sqnum.length - 1; i < maxLength; i++)
{
- cons2[i][gapIndex]++; // gap count
+ cons2GapCounts[i]++;
}
j++;
}
}
/**
- * Calculates the quality of the set of sequences
+ * Calculates the quality of the set of sequences over the given inclusive
+ * column range, using the specified substitution score matrix
*
- * @param startRes
- * Start residue
- * @param endRes
- * End residue
+ * @param startCol
+ * @param endCol
+ * @param scoreMatrix
*/
- public void findQuality(int startRes, int endRes)
+ protected void findQuality(int startCol, int endCol,
+ ScoreMatrix scoreMatrix)
{
- quality = new Vector<Double>();
+ quality = new Vector<>();
double max = -Double.MAX_VALUE;
- ScoreMatrix blosum62 = ScoreModels.getInstance().getBlosum62();
- float[][] blosumScores = blosum62.getMatrix();
+ float[][] scores = scoreMatrix.getMatrix();
- // Loop over columns // JBPNote Profiling info
- // long ts = System.currentTimeMillis();
- // long te = System.currentTimeMillis();
- percentIdentity(blosum62);
+ percentIdentity(scoreMatrix);
int size = seqNums.size();
int[] lengths = new int[size];
- double tot, bigtot, sr, tmp;
- double[] x, xx;
- int l, j, i, ii, i2, k, seqNum;
- for (l = 0; l < size; l++)
+ for (int l = 0; l < size; l++)
{
lengths[l] = seqNums.elementAt(l).length - 1;
}
- final int symbolCount = blosum62.getSize();
- int gapIndex = blosum62.getGapIndex();
+ final int symbolCount = scoreMatrix.getSize();
- for (j = startRes; j <= endRes; j++)
+ for (int j = startCol; j <= endCol; j++)
{
- bigtot = 0;
+ double bigtot = 0;
// First Xr = depends on column only
- x = new double[symbolCount];
+ double[] x = new double[symbolCount];
- for (ii = 0; ii < symbolCount; ii++)
+ for (int ii = 0; ii < symbolCount; ii++)
{
x[ii] = 0;
- for (i2 = 0; i2 < symbolCount; i2++)
+ /*
+ * todo JAL-728 currently assuming last symbol in matrix is * for gap
+ * (which we ignore as counted separately); true for BLOSUM62 but may
+ * not be once alternative matrices are supported
+ */
+ for (int i2 = 0; i2 < symbolCount - 1; i2++)
{
- x[ii] += (((double) cons2[j][i2] * blosumScores[ii][i2]) + 4);
+ x[ii] += (((double) cons2[j][i2] * scores[ii][i2]) + 4D);
}
+ x[ii] += 4D + cons2GapCounts[j] * scoreMatrix.getMinimumScore();
x[ii] /= size;
}
// Now calculate D for each position and sum
- for (k = 0; k < size; k++)
+ for (int k = 0; k < size; k++)
{
- tot = 0;
- xx = new double[symbolCount];
- seqNum = (j < lengths[k]) ? seqNums.elementAt(k)[j + 1] : gapIndex;
- // Sequence, or gap at the end
+ double tot = 0;
+ double[] xx = new double[symbolCount];
+ // sequence character index, or implied gap if sequence too short
+ int seqNum = (j < lengths[k]) ? seqNums.elementAt(k)[j + 1]
+ : GAP_INDEX;
- // This is a loop over r
- for (i = 0; i < symbolCount - 1; i++)
+ for (int i = 0; i < symbolCount - 1; i++)
{
- sr = 0;
-
- sr = (double) blosumScores[i][seqNum] + 4;
+ double sr = 4D;
+ if (seqNum == GAP_INDEX)
+ {
+ sr += scoreMatrix.getMinimumScore();
+ }
+ else
+ {
+ sr += scores[i][seqNum];
+ }
- // Calculate X with another loop over residues
- // System.out.println("Xi " + i + " " + x[i] + " " + sr);
xx[i] = x[i] - sr;
tot += (xx[i] * xx[i]);
bigtot += Math.sqrt(tot);
}
- // This is the quality for one column
- if (max < bigtot)
- {
- max = bigtot;
- }
+ max = Math.max(max, bigtot);
- // bigtot = bigtot * (size-cons2[j][23])/size;
- quality.addElement(new Double(bigtot));
-
- // Need to normalize by gaps
+ quality.addElement(Double.valueOf(bigtot));
}
double newmax = -Double.MAX_VALUE;
- for (j = startRes; j <= endRes; j++)
+ for (int j = startCol; j <= endCol; j++)
{
- tmp = quality.elementAt(j).doubleValue();
+ double tmp = quality.elementAt(j).doubleValue();
// tmp = ((max - tmp) * (size - cons2[j][23])) / size;
- tmp = ((max - tmp) * (size - cons2[j][gapIndex])) / size;
+ tmp = ((max - tmp) * (size - cons2GapCounts[j])) / size;
// System.out.println(tmp+ " " + j);
- quality.setElementAt(new Double(tmp), j);
+ quality.setElementAt(Double.valueOf(tmp), j);
if (tmp > newmax)
{
}
}
- // System.out.println("Quality " + s);
- qualityRange[0] = 0D;
- qualityRange[1] = newmax;
+ qualityMinimum = 0D;
+ qualityMaximum = newmax;
}
/**
* Complete the given consensus and quuality annotation rows. Note: currently
- * this method will enlarge the given annotation row if it is too small,
- * otherwise will leave its length unchanged.
+ * this method will reallocate the given annotation row if it is different to
+ * the calculated width, otherwise will leave its length unchanged.
*
* @param conservation
* conservation annotation row
public void completeAnnotations(AlignmentAnnotation conservation,
AlignmentAnnotation quality2, int istart, int alWidth)
{
- char[] sequence = getConsSequence().getSequence();
- float minR;
- float minG;
- float minB;
- float maxR;
- float maxG;
- float maxB;
- minR = 0.3f;
- minG = 0.0f;
- minB = 0f;
- maxR = 1.0f - minR;
- maxG = 0.9f - minG;
- maxB = 0f - minB; // scalable range for colouring both Conservation and
- // Quality
+ SequenceI cons = getConsSequence();
+
+ /*
+ * colour scale for Conservation and Quality;
+ */
+ float minR = 0.3f;
+ float minG = 0.0f;
+ float minB = 0f;
+ float maxR = 1.0f - minR;
+ float maxG = 0.9f - minG;
+ float maxB = 0f - minB;
float min = 0f;
float max = 11f;
float qmin = 0f;
float qmax = 0f;
- char c;
-
if (conservation != null && conservation.annotations != null
- && conservation.annotations.length < alWidth)
+ && conservation.annotations.length != alWidth)
{
conservation.annotations = new Annotation[alWidth];
}
if (quality2 != null)
{
- quality2.graphMax = (float) qualityRange[1];
+ quality2.graphMax = (float) qualityMaximum;
if (quality2.annotations != null
- && quality2.annotations.length < alWidth)
+ && quality2.annotations.length != alWidth)
{
quality2.annotations = new Annotation[alWidth];
}
- qmin = (float) qualityRange[0];
- qmax = (float) qualityRange[1];
+ qmin = (float) qualityMinimum;
+ qmax = (float) qualityMaximum;
}
for (int i = istart; i < alWidth; i++)
{
float value = 0;
- c = sequence[i];
+ char c = cons.getCharAt(i);
if (Character.isDigit(c))
{
float vprop = value - min;
vprop /= max;
int consp = i - start;
- String conssym = (value > 0 && consp > -1 && consp < consSymbs.length) ? consSymbs[consp]
- : "";
+ String conssym = (value > 0 && consp > -1
+ && consp < consSymbs.length) ? consSymbs[consp] : "";
conservation.annotations[i] = new Annotation(String.valueOf(c),
- conssym, ' ', value, new Color(minR + (maxR * vprop), minG
- + (maxG * vprop), minB + (maxB * vprop)));
+ conssym, ' ', value, new Color(minR + (maxR * vprop),
+ minG + (maxG * vprop), minB + (maxB * vprop)));
}
// Quality calc
value = quality.elementAt(i).floatValue();
float vprop = value - qmin;
vprop /= qmax;
- quality2.annotations[i] = new Annotation(" ",
- String.valueOf(value), ' ', value, new Color(minR
- + (maxR * vprop), minG + (maxG * vprop), minB
- + (maxB * vprop)));
+ String description = FORMAT_3DP.form(value);
+ quality2.annotations[i] = new Annotation(" ", description,
+ ' ', value, new Color(minR + (maxR * vprop),
+ minG + (maxG * vprop), minB + (maxB * vprop)));
}
}
}
*/
String getTooltip(int column)
{
- char[] sequence = getConsSequence().getSequence();
- char val = column < sequence.length ? sequence[column] : '-';
+ SequenceI cons = getConsSequence();
+ char val = column < cons.getLength() ? cons.getCharAt(column) : '-';
boolean hasConservation = val != '-' && val != '0';
int consp = column - start;
- String tip = (hasConservation && consp > -1 && consp < consSymbs.length) ? consSymbs[consp]
+ String tip = (hasConservation && consp > -1 && consp < consSymbs.length)
+ ? consSymbs[consp]
: "";
return tip;
}
git://source.jalview.org / jalview.git / blobdiff