X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=src%2Fjalview%2Fanalysis%2FStructureFrequency.java;h=c04df6c9ac56d71b71adcb1aacd73a36b92c7bd7;hb=4f77328104498504339216829abf5ea87e2791ec;hp=33a08e679bf2e3d5e69937b087d9bf347941077d;hpb=b57a02c25e335d033c97f8a6bacd6b54f62bd2b6;p=jalview.git diff --git a/src/jalview/analysis/StructureFrequency.java b/src/jalview/analysis/StructureFrequency.java index 33a08e6..c04df6c 100644 --- a/src/jalview/analysis/StructureFrequency.java +++ b/src/jalview/analysis/StructureFrequency.java @@ -1,26 +1,33 @@ /* - * 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. - * - * 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 + * 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. + * + * 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.Hashtable; /** * Takes in a vector or array of sequences and column start and column end and @@ -33,6 +40,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 +67,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; @@ -80,22 +93,27 @@ public class StructureFrequency * @param rnaStruc */ public static final void calculate(SequenceI[] sequences, int start, - int end, Hashtable[] result, boolean profile, + int end, Hashtable[] result, boolean profile, AlignmentAnnotation rnaStruc) { - Hashtable residueHash; + + 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 +131,7 @@ public class StructureFrequency s = '-'; } - if (s != '(') + if (!Rna.isOpeningParenthesis(s)) { if (s == '-') { @@ -123,63 +141,100 @@ public class StructureFrequency else { bpEnd = findPair(rna, i); + if (bpEnd > -1) { - for (j = 0; j < jSize; j++) // foreach row + for (int j = 0; j < jSize; j++) // foreach row { if (sequences[j] == null) { - System.err - .println("WARNING: Consensus skipping null sequence - possible race condition."); + System.err.println( + "WARNING: Consensus skipping null sequence - possible race condition."); continue; } - c = sequences[j].getCharAt(i); - { - // standard representation for gaps in sequence and structure - if (c == '.' || c == ' ') - { - c = '-'; - } + c = sequences[j].getCharAt(i); + cEnd = sequences[j].getCharAt(bpEnd); - if (c == '-') - { - values['-']++; - continue; - } - cEnd = sequences[j].getCharAt(bpEnd); - if (checkBpType(c, cEnd)) + if (Comparison.isGap(c) || Comparison.isGap(cEnd)) + { + values['-']++; + continue; + } + nongap++; + /* + * ensure upper-case for counting purposes + */ + if ('a' <= c && 'z' >= c) + { + c += 'A' - 'a'; + } + if ('a' <= cEnd && 'z' >= cEnd) + { + cEnd += 'A' - 'a'; + } + if (Rna.isCanonicalOrWobblePair(c, cEnd)) + { + canonicalOrWobblePairCount++; + if (Rna.isCanonicalPair(c, cEnd)) { - values['(']++; // H means it's a helix (structured) + canonical++; } - pairs[c][cEnd]++; - - maxResidue = "("; } + else + { + otherPairCount++; + } + pairs[c][cEnd]++; } } - // 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['(']; - - residueHash.put(MAXCOUNT, new Integer(count)); + 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, Integer.valueOf(count)); residueHash.put(MAXRESIDUE, maxResidue); percentage = ((float) count * 100) / jSize; - residueHash.put(PID_GAPS, new Float(percentage)); + residueHash.put(PID_GAPS, Float.valueOf(percentage)); + + percentage = ((float) count * 100) / nongap; + residueHash.put(PID_NOGAPS, Float.valueOf(percentage)); - // percentage = ((float) count * 100) / (float) nongap; - // residueHash.put(PID_NOGAPS, new Float(percentage)); if (result[i] == null) { result[i] = residueHash; @@ -187,100 +242,36 @@ 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 = ")"; - + residueHash = new Hashtable<>(); 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); } - residueHash.put(MAXCOUNT, new Integer(count)); + residueHash.put(MAXCOUNT, Integer.valueOf(count)); residueHash.put(MAXRESIDUE, maxResidue); percentage = ((float) count * 100) / jSize; - residueHash.put(PID_GAPS, new Float(percentage)); + residueHash.put(PID_GAPS, Float.valueOf(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, Float.valueOf(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; } /** @@ -296,9 +287,9 @@ public class StructureFrequency * @param includeAllConsSymbols */ public static void completeConsensus(AlignmentAnnotation consensus, - Hashtable[] hconsensus, int iStart, int width, + 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 +299,22 @@ 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++) { - Hashtable hci; + Hashtable hci; if (i >= hconsensus.length || ((hci = hconsensus[i]) == null)) { // happens if sequences calculated over were shorter than alignment @@ -361,15 +365,14 @@ 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 }; + ca[x] = new int[] { c, d }; vl[x] = pairs[c][d]; x++; } @@ -377,14 +380,18 @@ public class StructureFrequency 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 = (vl[c] * 100f / 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++; } @@ -394,7 +401,7 @@ public class StructureFrequency } else { - mouseOver += ((int) value + "%"); + mouseOver += (fmt.form(value) + "%"); } consensus.annotations[i] = new Annotation(maxRes, mouseOver, ' ', value); @@ -407,66 +414,57 @@ public class StructureFrequency * @param hconsensus * @return profile of the given column */ - public static int[] extractProfile(Hashtable hconsensus, + public static int[] extractProfile(Hashtable hconsensus, boolean ignoreGapsInConsensusCalculation) { - int[] rtnval = new int[74]; // 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 }; + 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) (vl[c] * 100f / 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; - } - - public static void main(String args[]) - { - // Short test to see if checkBpType works - ArrayList test = new ArrayList(); - test.add("A"); - test.add("c"); - test.add("g"); - test.add("T"); - test.add("U"); - for (String i : test) - { - for (String j : test) - { - System.out.println(i + "-" + j + ": " - + StructureFrequency.checkBpType(i.charAt(0), j.charAt(0))); - } - } + // 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; } }