From: jprocter Date: Wed, 11 Oct 2006 09:53:46 +0000 (+0000) Subject: fixes exception (array index out of bounds - due to off by one error) when applying... X-Git-Tag: Release_2_2~285 X-Git-Url: http://source.jalview.org/gitweb/?a=commitdiff_plain;h=b70a873b1f44dd388a1fda4c8cd9c30460e6679b;hp=7ae574fd34bd7ffd1ea23be7602d6b3dc57aa49f;p=jalview.git fixes exception (array index out of bounds - due to off by one error) when applying conservation colouring to a group --- diff --git a/src/jalview/analysis/Conservation.java b/src/jalview/analysis/Conservation.java index d1b508b..0af0689 100755 --- a/src/jalview/analysis/Conservation.java +++ b/src/jalview/analysis/Conservation.java @@ -1,539 +1,540 @@ -/* -* Jalview - A Sequence Alignment Editor and Viewer -* Copyright (C) 2006 AM Waterhouse, J Procter, G Barton, M Clamp, S Searle -* -* This program is free software; you can redistribute it and/or -* modify it under the terms of the GNU General Public License -* as published by the Free Software Foundation; either version 2 -* of the License, or (at your option) any later version. -* -* This program is distributed in the hope that it will be useful, -* but WITHOUT ANY WARRANTY; without even the implied warranty of -* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -* GNU General Public License for more details. -* -* You should have received a copy of the GNU General Public License -* along with this program; if not, write to the Free Software -* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA -*/ -package jalview.analysis; - -import jalview.datamodel.*; - -import java.util.*; - - -/** - * Calculates conservation values for a given set of sequences - * - * @author $author$ - * @version $Revision$ - */ -public class Conservation -{ - SequenceI [] sequences; - int start; - int end; - Vector seqNums; // vector of int vectors where first is sequence checksum - int maxLength = 0; // used by quality calcs - boolean seqNumsChanged = false; // updated after any change via calcSeqNum; - Hashtable [] total; - - /** Stores calculated quality values */ - public Vector quality; - - /** Stores maximum and minimum values of quality values */ - public Double[] qualityRange = new Double[2]; - String consString = ""; - Sequence consSequence; - Hashtable propHash; - int threshold; - String name = ""; - int[][] cons2; - - /** - * Creates a new Conservation object. - * - * @param name Name of conservation - * @param propHash DOCUMENT ME! - * @param threshold to count the residues in residueHash(). commonly used value is 3 - * @param sequences sequences to be used in calculation - * @param start start residue position - * @param end end residue position - */ - public Conservation(String name, Hashtable propHash, int threshold, - Vector sequences, int start, int end) - { - - this.name = name; - this.propHash = propHash; - this.threshold = threshold; - this.start = start; - this.end = end; - - - int s, sSize = sequences.size(); - SequenceI[] sarray = new SequenceI[sSize]; - this.sequences = sarray; - - for (s = 0; s < sSize; s++) - { - sarray[s] = (SequenceI) sequences.elementAt(s); - if(sarray[s].getLength()>maxLength) - maxLength = sarray[s].getLength(); - } - } - - - /** - * DOCUMENT ME! - * - * @param i DOCUMENT ME! - */ - private void calcSeqNum(int i) - { - String sq = null; // for dumb jbuilder not-inited exception warning - int[] sqnum = null; - - int sSize = sequences.length; - - if ((i > -1) && (i < sSize)) - { - sq = sequences[i].getSequence(); - - if (seqNums.size() <= i) - { - seqNums.addElement(new int[sq.length() + 1]); - } - - if (sq.hashCode() != ((int[]) seqNums.elementAt(i))[0]) - { - int j; - int len; - seqNumsChanged = true; - len = sq.length(); - - if (maxLength < len) - { - maxLength = len; - } - - sqnum = new int[len + 1]; // better to always make a new array - sequence can change its length - sqnum[0] = sq.hashCode(); - - for (j = 1; j <= len; j++) - { - sqnum[j] = jalview.schemes.ResidueProperties.aaIndex[sq.charAt(j-1)]; - } - - - seqNums.setElementAt(sqnum, i); - } - else - System.out.println("NEVER THE EXCEPTION"); - } - else - { - // JBPNote INFO level debug - System.err.println( - "ERROR: calcSeqNum called with out of range sequence index for Alignment\n"); - } - } - - /** - * Calculates the conservation values for given set of sequences - */ - public void calculate() - { - Hashtable resultHash, ht; - int thresh, j, jSize = sequences.length; - int[] values; // Replaces residueHash - String type, res=null; - char c; - Enumeration enumeration2; - - total = new Hashtable[maxLength]; - - for (int i = start; i <= end; i++) - { - values = new int[132]; - - for (j = 0; j < jSize; j++) - { - if (sequences[j].getLength() > i) - { - c = sequences[j].getCharAt(i); - - // No need to check if its a '-' - if (c == '.' || c == ' ') - c = '-'; - - if ('a' <= c && c <= 'z') - { - c -= (32);// 32 = 'a' - 'A' - } - - values[c]++; - } - else - { - values['-']++; - } - } - - //What is the count threshold to count the residues in residueHash() - thresh = (threshold * (jSize)) / 100; - - //loop over all the found residues - resultHash = new Hashtable(); - for (int v = '-'; v < 'Z'; v++) - { - - if (values[v] > thresh) - { - res = String.valueOf( (char) v); - - //Now loop over the properties - enumeration2 = propHash.keys(); - - while (enumeration2.hasMoreElements()) - { - type = (String) enumeration2.nextElement(); - ht = (Hashtable) propHash.get(type); - - //Have we ticked this before? - if (!resultHash.containsKey(type)) - { - if (ht.containsKey(res)) - { - resultHash.put(type, ht.get(res)); - } - else - { - resultHash.put(type, ht.get("-")); - } - } - else if (((Integer) resultHash.get(type)).equals( - (Integer) ht.get(res)) == false) - { - resultHash.put(type, new Integer(-1)); - } - } - } - } - - total[i] = resultHash; - } - } - - - /*** - * countConsNGaps - * returns gap count in int[0], and conserved residue count in int[1] - */ - public int[] countConsNGaps(int j) - { - int count = 0; - int cons = 0; - int nres = 0; - int[] r = new int[2]; - char f = '$'; - int i, iSize = sequences.length; - char c; - - for (i = 0; i < iSize; i++) - { - if (j >= sequences[i].getLength()) - { - count++; - continue; - } - - c = sequences[i].getCharAt(j); // gaps do not have upper/lower case - - if (jalview.util.Comparison.isGap((c))) - { - count++; - } - else - { - nres++; - - if (nres == 1) - { - f = c; - cons++; - } - else if (f == c) - { - cons++; - } - } - } - - r[0] = (nres == cons) ? 1 : 0; - r[1] = count; - - return r; - } - - /** - * Calculates the conservation sequence - * - * @param consflag if true, poitiveve conservation; false calculates negative conservation - * @param percentageGaps commonly used value is 25 - */ - public void verdict(boolean consflag, float percentageGaps) - { - StringBuffer consString = new StringBuffer(); - String type; - Integer result; - int[] gapcons; - int totGaps, count; - float pgaps; - Hashtable resultHash ; - Enumeration enumeration; - - - for (int i = start; i <= end; i++) - { - gapcons = countConsNGaps(i); - totGaps = gapcons[1]; - pgaps = ((float) totGaps * 100) / (float) sequences.length; - - if (percentageGaps > pgaps) - { - resultHash = total[i - start]; - - //Now find the verdict - count = 0; - enumeration = resultHash.keys(); - - while (enumeration.hasMoreElements()) - { - type = (String) enumeration.nextElement(); - result = (Integer) resultHash.get(type); - - //Do we want to count +ve conservation or +ve and -ve cons.? - if (consflag) - { - if (result.intValue() == 1) - { - count++; - } - } - else - { - if (result.intValue() != -1) - { - count++; - } - } - } - - if (count < 10) - { - consString.append(count); // Conserved props!=Identity - } - else - { - consString.append((gapcons[0] == 1) ? "*" : "+"); - } - } - else - { - consString.append("-"); - } - } - - consSequence = new Sequence(name, consString.toString(), start, end); - } - - /** - * - * - * @return Conservation sequence - */ - public Sequence getConsSequence() - { - return consSequence; - } - - // From Alignment.java in jalview118 - public void findQuality() - { - findQuality(0, maxLength - 1); - } - - /** - * DOCUMENT ME! - */ - private void percentIdentity2() - { - seqNums = new Vector(); - // calcSeqNum(s); - int i = 0, iSize = sequences.length; - //Do we need to calculate this again? - for (i = 0; i < iSize; i++) - { - calcSeqNum(i); - } - - - if ((cons2 == null) || seqNumsChanged) - { - cons2 = new int[maxLength][24]; - - // Initialize the array - for (int j = 0; j < 24; j++) - { - for (i = 0; i < maxLength; i++) - { - cons2[i][j] = 0; - } - } - - int[] sqnum; - int j = 0; - - while (j < sequences.length) - { - sqnum = (int[]) seqNums.elementAt(j); - - for (i = 1; i < sqnum.length; i++) - { - cons2[i - 1][sqnum[i]]++; - } - - for (i = sqnum.length - 1; i < maxLength; i++) - { - cons2[i][23]++; // gap count - } - - j++; - } - - // unnecessary ? - - /* for (int i=start; i <= end; i++) { - int max = -1000; - int maxi = -1; - int maxj = -1; - - for (int j=0;j<24;j++) { - if (cons2[i][j] > max) { - max = cons2[i][j]; - maxi = i; - maxj = j; - } - - } - } */ - } - } - - /** - * Calculates the quality of the set of sequences - * - * @param start Start residue - * @param end End residue - */ - public void findQuality(int start, int end) - { - quality = new Vector(); - - double max = -10000; - int[][] BLOSUM62 = jalview.schemes.ResidueProperties.getBLOSUM62(); - - //Loop over columns // JBPNote Profiling info - //long ts = System.currentTimeMillis(); - //long te = System.currentTimeMillis(); - percentIdentity2(); - - 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++) - lengths[l] = ((int[]) seqNums.elementAt(l)).length - 1; - - for (j = start; j <= end; j++) - { - bigtot = 0; - - // First Xr = depends on column only - x = new double[24]; - - for (ii = 0; ii < 24; ii++) - { - x[ii] = 0; - - for (i2 = 0; i2 < 24; i2++) - { - x[ii] += ( ( (double) cons2[j][i2] * BLOSUM62[ii][i2]) + - 4); - } - - x[ii] /= size; - } - - // Now calculate D for each position and sum - for (k = 0; k < size; k++) - { - tot = 0; - xx = new double[24]; - seqNum = (j < lengths[k]) - ? ((int[]) seqNums.elementAt(k))[j + 1] : 23; // Sequence, or gap at the end - - // This is a loop over r - for (i = 0; i < 23; i++) - { - sr = 0; - - sr = (double) BLOSUM62[i][seqNum] + 4; - - //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; - } - - // bigtot = bigtot * (size-cons2[j][23])/size; - quality.addElement(new Double(bigtot)); - - // Need to normalize by gaps - } - - double newmax = -10000; - - for (j = start; j <= end; j++) - { - tmp = ((Double) quality.elementAt(j)).doubleValue(); - tmp = ((max - tmp) * (size - cons2[j][23])) / size; - - // System.out.println(tmp+ " " + j); - quality.setElementAt(new Double(tmp), j); - - if (tmp > newmax) - { - newmax = tmp; - } - } - - // System.out.println("Quality " + s); - qualityRange[0] = new Double(0); - qualityRange[1] = new Double(newmax); - } -} +/* +* Jalview - A Sequence Alignment Editor and Viewer +* Copyright (C) 2006 AM Waterhouse, J Procter, G Barton, M Clamp, S Searle +* +* This program is free software; you can redistribute it and/or +* modify it under the terms of the GNU General Public License +* as published by the Free Software Foundation; either version 2 +* of the License, or (at your option) any later version. +* +* This program is distributed in the hope that it will be useful, +* but WITHOUT ANY WARRANTY; without even the implied warranty of +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +* GNU General Public License for more details. +* +* You should have received a copy of the GNU General Public License +* along with this program; if not, write to the Free Software +* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA +*/ +package jalview.analysis; + +import jalview.datamodel.*; + +import java.util.*; + + +/** + * Calculates conservation values for a given set of sequences + * + * @author $author$ + * @version $Revision$ + */ +public class Conservation +{ + SequenceI [] sequences; + int start; + int end; + Vector seqNums; // vector of int vectors where first is sequence checksum + int maxLength = 0; // used by quality calcs + boolean seqNumsChanged = false; // updated after any change via calcSeqNum; + Hashtable [] total; + + /** Stores calculated quality values */ + public Vector quality; + + /** Stores maximum and minimum values of quality values */ + public Double[] qualityRange = new Double[2]; + String consString = ""; + Sequence consSequence; + Hashtable propHash; + int threshold; + String name = ""; + int[][] cons2; + + /** + * Creates a new Conservation object. + * + * @param name Name of conservation + * @param propHash DOCUMENT ME! + * @param threshold to count the residues in residueHash(). commonly used value is 3 + * @param sequences sequences to be used in calculation + * @param start start residue position + * @param end end residue position + */ + public Conservation(String name, Hashtable propHash, int threshold, + Vector sequences, int start, int end) + { + + this.name = name; + this.propHash = propHash; + this.threshold = threshold; + this.start = start; + this.end = end; + + maxLength=end-start+1; // default width includes bounds of calculation + + int s, sSize = sequences.size(); + SequenceI[] sarray = new SequenceI[sSize]; + this.sequences = sarray; + + for (s = 0; s < sSize; s++) + { + sarray[s] = (SequenceI) sequences.elementAt(s); + if(sarray[s].getLength()>maxLength) + maxLength = sarray[s].getLength(); + } + } + + + /** + * DOCUMENT ME! + * + * @param i DOCUMENT ME! + */ + private void calcSeqNum(int i) + { + String sq = null; // for dumb jbuilder not-inited exception warning + int[] sqnum = null; + + int sSize = sequences.length; + + if ((i > -1) && (i < sSize)) + { + sq = sequences[i].getSequence(); + + if (seqNums.size() <= i) + { + seqNums.addElement(new int[sq.length() + 1]); + } + + if (sq.hashCode() != ((int[]) seqNums.elementAt(i))[0]) + { + int j; + int len; + seqNumsChanged = true; + len = sq.length(); + + if (maxLength < len) + { + maxLength = len; + } + + sqnum = new int[len + 1]; // better to always make a new array - sequence can change its length + sqnum[0] = sq.hashCode(); + + for (j = 1; j <= len; j++) + { + sqnum[j] = jalview.schemes.ResidueProperties.aaIndex[sq.charAt(j-1)]; + } + + + seqNums.setElementAt(sqnum, i); + } + else + System.out.println("NEVER THE EXCEPTION"); + } + else + { + // JBPNote INFO level debug + System.err.println( + "ERROR: calcSeqNum called with out of range sequence index for Alignment\n"); + } + } + + /** + * Calculates the conservation values for given set of sequences + */ + public void calculate() + { + Hashtable resultHash, ht; + int thresh, j, jSize = sequences.length; + int[] values; // Replaces residueHash + String type, res=null; + char c; + Enumeration enumeration2; + + total = new Hashtable[maxLength]; + + for (int i = start; i <= end; i++) + { + values = new int[132]; + + for (j = 0; j < jSize; j++) + { + if (sequences[j].getLength() > i) + { + c = sequences[j].getCharAt(i); + + // No need to check if its a '-' + if (c == '.' || c == ' ') + c = '-'; + + if ('a' <= c && c <= 'z') + { + c -= (32);// 32 = 'a' - 'A' + } + + values[c]++; + } + else + { + values['-']++; + } + } + + //What is the count threshold to count the residues in residueHash() + thresh = (threshold * (jSize)) / 100; + + //loop over all the found residues + resultHash = new Hashtable(); + for (int v = '-'; v < 'Z'; v++) + { + + if (values[v] > thresh) + { + res = String.valueOf( (char) v); + + //Now loop over the properties + enumeration2 = propHash.keys(); + + while (enumeration2.hasMoreElements()) + { + type = (String) enumeration2.nextElement(); + ht = (Hashtable) propHash.get(type); + + //Have we ticked this before? + if (!resultHash.containsKey(type)) + { + if (ht.containsKey(res)) + { + resultHash.put(type, ht.get(res)); + } + else + { + resultHash.put(type, ht.get("-")); + } + } + else if (((Integer) resultHash.get(type)).equals( + (Integer) ht.get(res)) == false) + { + resultHash.put(type, new Integer(-1)); + } + } + } + } + + total[i-start] = resultHash; + } + } + + + /*** + * countConsNGaps + * returns gap count in int[0], and conserved residue count in int[1] + */ + public int[] countConsNGaps(int j) + { + int count = 0; + int cons = 0; + int nres = 0; + int[] r = new int[2]; + char f = '$'; + int i, iSize = sequences.length; + char c; + + for (i = 0; i < iSize; i++) + { + if (j >= sequences[i].getLength()) + { + count++; + continue; + } + + c = sequences[i].getCharAt(j); // gaps do not have upper/lower case + + if (jalview.util.Comparison.isGap((c))) + { + count++; + } + else + { + nres++; + + if (nres == 1) + { + f = c; + cons++; + } + else if (f == c) + { + cons++; + } + } + } + + r[0] = (nres == cons) ? 1 : 0; + r[1] = count; + + return r; + } + + /** + * Calculates the conservation sequence + * + * @param consflag if true, poitiveve conservation; false calculates negative conservation + * @param percentageGaps commonly used value is 25 + */ + public void verdict(boolean consflag, float percentageGaps) + { + StringBuffer consString = new StringBuffer(); + String type; + Integer result; + int[] gapcons; + int totGaps, count; + float pgaps; + Hashtable resultHash ; + Enumeration enumeration; + + + for (int i = start; i <= end; i++) + { + gapcons = countConsNGaps(i); + totGaps = gapcons[1]; + pgaps = ((float) totGaps * 100) / (float) sequences.length; + + if (percentageGaps > pgaps) + { + resultHash = total[i - start]; + + //Now find the verdict + count = 0; + enumeration = resultHash.keys(); + + while (enumeration.hasMoreElements()) + { + type = (String) enumeration.nextElement(); + result = (Integer) resultHash.get(type); + + //Do we want to count +ve conservation or +ve and -ve cons.? + if (consflag) + { + if (result.intValue() == 1) + { + count++; + } + } + else + { + if (result.intValue() != -1) + { + count++; + } + } + } + + if (count < 10) + { + consString.append(count); // Conserved props!=Identity + } + else + { + consString.append((gapcons[0] == 1) ? "*" : "+"); + } + } + else + { + consString.append("-"); + } + } + + consSequence = new Sequence(name, consString.toString(), start, end); + } + + /** + * + * + * @return Conservation sequence + */ + public Sequence getConsSequence() + { + return consSequence; + } + + // From Alignment.java in jalview118 + public void findQuality() + { + findQuality(0, maxLength - 1); + } + + /** + * DOCUMENT ME! + */ + private void percentIdentity2() + { + seqNums = new Vector(); + // calcSeqNum(s); + int i = 0, iSize = sequences.length; + //Do we need to calculate this again? + for (i = 0; i < iSize; i++) + { + calcSeqNum(i); + } + + + if ((cons2 == null) || seqNumsChanged) + { + cons2 = new int[maxLength][24]; + + // Initialize the array + for (int j = 0; j < 24; j++) + { + for (i = 0; i < maxLength; i++) + { + cons2[i][j] = 0; + } + } + + int[] sqnum; + int j = 0; + + while (j < sequences.length) + { + sqnum = (int[]) seqNums.elementAt(j); + + for (i = 1; i < sqnum.length; i++) + { + cons2[i - 1][sqnum[i]]++; + } + + for (i = sqnum.length - 1; i < maxLength; i++) + { + cons2[i][23]++; // gap count + } + + j++; + } + + // unnecessary ? + + /* for (int i=start; i <= end; i++) { + int max = -1000; + int maxi = -1; + int maxj = -1; + + for (int j=0;j<24;j++) { + if (cons2[i][j] > max) { + max = cons2[i][j]; + maxi = i; + maxj = j; + } + + } + } */ + } + } + + /** + * Calculates the quality of the set of sequences + * + * @param start Start residue + * @param end End residue + */ + public void findQuality(int start, int end) + { + quality = new Vector(); + + double max = -10000; + int[][] BLOSUM62 = jalview.schemes.ResidueProperties.getBLOSUM62(); + + //Loop over columns // JBPNote Profiling info + //long ts = System.currentTimeMillis(); + //long te = System.currentTimeMillis(); + percentIdentity2(); + + 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++) + lengths[l] = ((int[]) seqNums.elementAt(l)).length - 1; + + for (j = start; j <= end; j++) + { + bigtot = 0; + + // First Xr = depends on column only + x = new double[24]; + + for (ii = 0; ii < 24; ii++) + { + x[ii] = 0; + + for (i2 = 0; i2 < 24; i2++) + { + x[ii] += ( ( (double) cons2[j][i2] * BLOSUM62[ii][i2]) + + 4); + } + + x[ii] /= size; + } + + // Now calculate D for each position and sum + for (k = 0; k < size; k++) + { + tot = 0; + xx = new double[24]; + seqNum = (j < lengths[k]) + ? ((int[]) seqNums.elementAt(k))[j + 1] : 23; // Sequence, or gap at the end + + // This is a loop over r + for (i = 0; i < 23; i++) + { + sr = 0; + + sr = (double) BLOSUM62[i][seqNum] + 4; + + //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; + } + + // bigtot = bigtot * (size-cons2[j][23])/size; + quality.addElement(new Double(bigtot)); + + // Need to normalize by gaps + } + + double newmax = -10000; + + for (j = start; j <= end; j++) + { + tmp = ((Double) quality.elementAt(j)).doubleValue(); + tmp = ((max - tmp) * (size - cons2[j][23])) / size; + + // System.out.println(tmp+ " " + j); + quality.setElementAt(new Double(tmp), j); + + if (tmp > newmax) + { + newmax = tmp; + } + } + + // System.out.println("Quality " + s); + qualityRange[0] = new Double(0); + qualityRange[1] = new Double(newmax); + } +}