X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=src%2Fjalview%2Fanalysis%2FStructureFrequency.java;h=7b0858d502d70a7a5aac3582478afcc6cea7f1e0;hb=37de9310bec3501cbc6381e0c3dcb282fcaad812;hp=aa1e2774beed6ad0f0ac5d3ef36c2758fd8e76a5;hpb=bbfdb203ce3a889600caa52478478b75ab03809f;p=jalview.git diff --git a/src/jalview/analysis/StructureFrequency.java b/src/jalview/analysis/StructureFrequency.java index aa1e277..7b0858d 100644 --- a/src/jalview/analysis/StructureFrequency.java +++ b/src/jalview/analysis/StructureFrequency.java @@ -1,26 +1,34 @@ /* - * Jalview - A Sequence Alignment Editor and Viewer (Version 2.7) - * Copyright (C) 2011 J Procter, AM Waterhouse, J Engelhardt, LM Lui, G Barton, M Clamp, S Searle + * Jalview - A Sequence Alignment Editor and Viewer ($$Version-Rel$$) + * Copyright (C) $$Year-Rel$$ The Jalview Authors * * This file is part of Jalview. * * Jalview 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 3 of the License, or (at your option) any later version. - * + * as published by the Free Software Foundation, either version 3 + * of the License, or (at your option) any later version. + * * Jalview 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 Jalview. If not, see . + * You should have received a copy of the GNU General Public License + * along with Jalview. If not, see . + * The Jalview Authors are detailed in the 'AUTHORS' file. */ - package jalview.analysis; -import java.util.*; +import jalview.datamodel.AlignmentAnnotation; +import jalview.datamodel.Annotation; +import jalview.datamodel.SequenceFeature; +import jalview.datamodel.SequenceI; +import jalview.util.Comparison; +import jalview.util.Format; -import jalview.datamodel.*; +import java.util.ArrayList; +import java.util.Hashtable; /** * Takes in a vector or array of sequences and column start and column end and @@ -33,6 +41,8 @@ import jalview.datamodel.*; */ public class StructureFrequency { + public static final int STRUCTURE_PROFILE_LENGTH = 74; + // No need to store 1000s of strings which are not // visible to the user. public static final String MAXCOUNT = "C"; @@ -58,11 +68,15 @@ public class StructureFrequency */ public static int findPair(SequenceFeature[] pairs, int indice) { + for (int i = 0; i < pairs.length; i++) { if (pairs[i].getBegin() == indice) + { + return pairs[i].getEnd(); + } } return -1; @@ -83,19 +97,24 @@ public class StructureFrequency int end, Hashtable[] result, boolean profile, AlignmentAnnotation rnaStruc) { + Hashtable residueHash; String maxResidue; - char[] seq, struc = rnaStruc.getRNAStruc().toCharArray(); + char[] struc = rnaStruc.getRNAStruc().toCharArray(); + SequenceFeature[] rna = rnaStruc._rnasecstr; char c, s, cEnd; - int count, nonGap = 0, i, bpEnd = -1, j, jSize = sequences.length; + int bpEnd = -1; + int jSize = sequences.length; int[] values; int[][] pairs; float percentage; - for (i = start; i < end; i++) // foreach column + for (int i = start; i < end; i++) // foreach column { - residueHash = new Hashtable(); + int canonicalOrWobblePairCount = 0, canonical = 0; + int otherPairCount = 0; + int nongap = 0; maxResidue = "-"; values = new int[255]; pairs = new int[255][255]; @@ -113,7 +132,7 @@ public class StructureFrequency s = '-'; } - if (s != '(') + if (!Rna.isOpeningParenthesis(s)) { if (s == '-') { @@ -122,64 +141,99 @@ public class StructureFrequency } else { - for (j = 0; j < jSize; j++) // foreach row - { - if (sequences[j] == null) - { - System.err - .println("WARNING: Consensus skipping null sequence - possible race condition."); - continue; - } - seq = sequences[j].getSequence(); + bpEnd = findPair(rna, i); - if (seq.length > i) + if (bpEnd > -1) + { + for (int j = 0; j < jSize; j++) // foreach row { - c = seq[i]; - - // standard representation for gaps in sequence and structure - if (c == '.' || c == ' ') + if (sequences[j] == null) { - c = '-'; + System.err + .println("WARNING: Consensus skipping null sequence - possible race condition."); + continue; } - if (c == '-') + c = sequences[j].getCharAt(i); + cEnd = sequences[j].getCharAt(bpEnd); + + if (Comparison.isGap(c) || Comparison.isGap(cEnd)) { values['-']++; continue; } - bpEnd = findPair(rna, i); - cEnd = seq[bpEnd]; - if (checkBpType(c, cEnd)) + nongap++; + /* + * ensure upper-case for counting purposes + */ + if ('a' <= c && 'z' >= c) + { + c += 'A' - 'a'; + } + if ('a' <= cEnd && 'z' >= cEnd) { - values['(']++; // H means it's a helix (structured) + cEnd += 'A' - 'a'; + } + if (Rna.isCanonicalOrWobblePair(c, cEnd)) + { + canonicalOrWobblePairCount++; + if (Rna.isCanonicalPair(c, cEnd)) + { + canonical++; + } + } + else + { + otherPairCount++; } pairs[c][cEnd]++; - - maxResidue = "("; } } - // nonGap++; } - // UPDATE this for new values + + residueHash = new Hashtable(); if (profile) { - residueHash.put(PROFILE, new int[][] - { values, new int[] - { jSize, (jSize - values['-']) } }); + // TODO 1-dim array with jsize in [0], nongapped in [1]; or Pojo + residueHash.put(PROFILE, new int[][] { values, + new int[] { jSize, (jSize - values['-']) } }); residueHash.put(PAIRPROFILE, pairs); } - - count = values['(']; - + values['('] = canonicalOrWobblePairCount; + values['['] = canonical; + values['{'] = otherPairCount; + /* + * the count is the number of valid pairs (as a percentage, determines + * the relative size of the profile logo) + */ + int count = canonicalOrWobblePairCount; + + /* + * display '(' if most pairs are canonical, or as + * '[' if there are more wobble pairs. + */ + if (canonicalOrWobblePairCount > 0 || otherPairCount > 0) + { + if (canonicalOrWobblePairCount >= otherPairCount) + { + maxResidue = (canonicalOrWobblePairCount - canonical) < canonical ? "(" + : "["; + } + else + { + maxResidue = "{"; + } + } residueHash.put(MAXCOUNT, new Integer(count)); residueHash.put(MAXRESIDUE, maxResidue); - percentage = ((float) count * 100) / (float) jSize; + percentage = ((float) count * 100) / jSize; residueHash.put(PID_GAPS, new Float(percentage)); - // percentage = ((float) count * 100) / (float) nongap; - // residueHash.put(PID_NOGAPS, new Float(percentage)); + percentage = ((float) count * 100) / nongap; + residueHash.put(PID_NOGAPS, new Float(percentage)); + if (result[i] == null) { result[i] = residueHash; @@ -187,16 +241,19 @@ public class StructureFrequency if (bpEnd > 0) { values[')'] = values['(']; + values[']'] = values['[']; + values['}'] = values['{']; values['('] = 0; + values['['] = 0; + values['{'] = 0; + maxResidue = maxResidue.equals("(") ? ")" + : maxResidue.equals("[") ? "]" : "}"; residueHash = new Hashtable(); - maxResidue = ")"; - if (profile) { - residueHash.put(PROFILE, new int[][] - { values, new int[] - { jSize, (jSize - values['-']) } }); + residueHash.put(PROFILE, new int[][] { values, + new int[] { jSize, (jSize - values['-']) } }); residueHash.put(PAIRPROFILE, pairs); } @@ -204,83 +261,15 @@ public class StructureFrequency residueHash.put(MAXCOUNT, new Integer(count)); residueHash.put(MAXRESIDUE, maxResidue); - percentage = ((float) count * 100) / (float) jSize; + percentage = ((float) count * 100) / jSize; residueHash.put(PID_GAPS, new Float(percentage)); - result[bpEnd] = residueHash; - } - } - } - - /** - * Method to check if a base-pair is a canonical or a wobble bp - * - * @param up - * 5' base - * @param down - * 3' base - * @return True if it is a canonical/wobble bp - */ - public static boolean checkBpType(char up, char down) - { - if (up > 'Z') - { - up -= 32; - } - if (down > 'Z') - { - down -= 32; - } + percentage = ((float) count * 100) / nongap; + residueHash.put(PID_NOGAPS, new Float(percentage)); - switch (up) - { - case 'A': - switch (down) - { - case 'T': - return true; - case 'U': - return true; - } - break; - case 'C': - switch (down) - { - case 'G': - return true; - } - break; - case 'T': - switch (down) - { - case 'A': - return true; - case 'G': - return true; - } - break; - case 'G': - switch (down) - { - case 'C': - return true; - case 'T': - return true; - case 'U': - return true; - } - break; - case 'U': - switch (down) - { - case 'A': - return true; - case 'G': - return true; + result[bpEnd] = residueHash; } - break; } - return false; } /** @@ -298,7 +287,7 @@ public class StructureFrequency public static void completeConsensus(AlignmentAnnotation consensus, Hashtable[] hconsensus, int iStart, int width, boolean ignoreGapsInConsensusCalculation, - boolean includeAllConsSymbols) + boolean includeAllConsSymbols, long nseq) { float tval, value; if (consensus == null || consensus.annotations == null @@ -308,9 +297,23 @@ public class StructureFrequency // initialised properly return; } + String fmtstr = "%3.1f"; + int precision = 2; + while (nseq > 100) + { + precision++; + nseq /= 10; + } + if (precision > 2) + { + fmtstr = "%" + (2 + precision) + "." + precision + "f"; + } + Format fmt = new Format(fmtstr); + for (int i = iStart; i < width; i++) { - if (i >= hconsensus.length) + Hashtable hci; + if (i >= hconsensus.length || ((hci = hconsensus[i]) == null)) { // happens if sequences calculated over were shorter than alignment // width @@ -318,30 +321,31 @@ public class StructureFrequency continue; } value = 0; + Float fv; if (ignoreGapsInConsensusCalculation) { - value = ((Float) hconsensus[i].get(StructureFrequency.PID_NOGAPS)) - .floatValue(); + fv = (Float) hci.get(StructureFrequency.PID_NOGAPS); } else { - value = ((Float) hconsensus[i].get(StructureFrequency.PID_GAPS)) - .floatValue(); + fv = (Float) hci.get(StructureFrequency.PID_GAPS); } - - String maxRes = hconsensus[i].get(StructureFrequency.MAXRESIDUE) - .toString(); - String mouseOver = hconsensus[i].get(StructureFrequency.MAXRESIDUE) - + " "; + if (fv == null) + { + consensus.annotations[i] = null; + // data has changed below us .. give up and + continue; + } + value = fv.floatValue(); + String maxRes = hci.get(StructureFrequency.MAXRESIDUE).toString(); + String mouseOver = hci.get(StructureFrequency.MAXRESIDUE) + " "; if (maxRes.length() > 1) { mouseOver = "[" + maxRes + "] "; maxRes = "+"; } - int[][] profile = (int[][]) hconsensus[i] - .get(StructureFrequency.PROFILE); - int[][] pairs = (int[][]) hconsensus[i] - .get(StructureFrequency.PAIRPROFILE); + int[][] profile = (int[][]) hci.get(StructureFrequency.PROFILE); + int[][] pairs = (int[][]) hci.get(StructureFrequency.PAIRPROFILE); if (pairs != null && includeAllConsSymbols) // Just responsible for the // tooltip @@ -349,7 +353,9 @@ public class StructureFrequency { mouseOver = ""; - /* TODO It's not sure what is the purpose of the alphabet and wheter it is useful for structure? + /* + * TODO It's not sure what is the purpose of the alphabet and wheter it + * is useful for structure? * * if (alphabet != null) { for (int c = 0; c < alphabet.length; c++) { * tval = ((float) profile[0][alphabet[c]]) 100f / (float) @@ -357,30 +363,33 @@ public class StructureFrequency * ((c == 0) ? "" : "; ") + alphabet[c] + " " + ((int) tval) + "%"; } } * else { */ - Object[] ca = new Object[625]; + int[][] ca = new int[625][]; float[] vl = new float[625]; int x = 0; for (int c = 65; c < 90; c++) { for (int d = 65; d < 90; d++) { - ca[x] = new int[] - { c, d }; - vl[x] = (float) pairs[c][d]; + ca[x] = new int[] { c, d }; + vl[x] = pairs[c][d]; x++; } } jalview.util.QuickSort.sort(vl, ca); int p = 0; + /* + * profile[1] is {total, ungappedTotal} + */ + final int divisor = profile[1][ignoreGapsInConsensusCalculation ? 1 + : 0]; for (int c = 624; c > 0; c--) { if (vl[c] > 0) { - tval = ((float) vl[c] * 100f / (float) profile[1][ignoreGapsInConsensusCalculation ? 1 - : 0]); - mouseOver += ((p == 0) ? "" : "; ") + (char) ((int[]) ca[c])[0] - + (char) ((int[]) ca[c])[1] + " " + ((int) tval) + "%"; + tval = (vl[c] * 100f / divisor); + mouseOver += ((p == 0) ? "" : "; ") + (char) ca[c][0] + + (char) ca[c][1] + " " + fmt.form(tval) + "%"; p++; } @@ -390,7 +399,7 @@ public class StructureFrequency } else { - mouseOver += ((int) value + "%"); + mouseOver += (fmt.form(value) + "%"); } consensus.annotations[i] = new Annotation(maxRes, mouseOver, ' ', value); @@ -406,45 +415,54 @@ public class StructureFrequency public static int[] extractProfile(Hashtable hconsensus, boolean ignoreGapsInConsensusCalculation) { - int[] rtnval = new int[52]; // 2*(5*5)+2 + int[] rtnval = new int[STRUCTURE_PROFILE_LENGTH]; // 2*(5*5)+2 int[][] profile = (int[][]) hconsensus.get(StructureFrequency.PROFILE); int[][] pairs = (int[][]) hconsensus .get(StructureFrequency.PAIRPROFILE); if (profile == null) + { return null; + } // TODO fix the object length, also do it in completeConsensus - Object[] ca = new Object[625]; + // Object[] ca = new Object[625]; + int[][] ca = new int[625][]; float[] vl = new float[625]; int x = 0; for (int c = 65; c < 90; c++) { for (int d = 65; d < 90; d++) { - ca[x] = new int[] - { c, d }; - vl[x] = (float) pairs[c][d]; + ca[x] = new int[] { c, d }; + vl[x] = pairs[c][d]; x++; } } jalview.util.QuickSort.sort(vl, ca); - rtnval[0] = 2; + int valuesCount = 0; rtnval[1] = 0; + int offset = 2; + final int divisor = profile[1][ignoreGapsInConsensusCalculation ? 1 : 0]; for (int c = 624; c > 0; c--) { if (vl[c] > 0) { - rtnval[rtnval[0]++] = ((int[]) ca[c])[0]; - rtnval[rtnval[0]++] = ((int[]) ca[c])[1]; - rtnval[rtnval[0]] = (int) ((float) vl[c] * 100f / (float) profile[1][ignoreGapsInConsensusCalculation ? 1 - : 0]); - rtnval[1]+=rtnval[rtnval[0]++]; + rtnval[offset++] = ca[c][0]; + rtnval[offset++] = ca[c][1]; + rtnval[offset] = (int) (vl[c] * 100f / divisor); + rtnval[1] += rtnval[offset++]; + valuesCount++; } } + rtnval[0] = valuesCount; - return rtnval; + // insert profile type code in position 0 + int[] result = new int[rtnval.length + 1]; + result[0] = AlignmentAnnotation.STRUCTURE_PROFILE; + System.arraycopy(rtnval, 0, result, 1, rtnval.length); + return result; } public static void main(String args[]) @@ -461,7 +479,7 @@ public class StructureFrequency for (String j : test) { System.out.println(i + "-" + j + ": " - + StructureFrequency.checkBpType(i.charAt(0), j.charAt(0))); + + Rna.isCanonicalOrWobblePair(i.charAt(0), j.charAt(0))); } } }