2 * Jalview - A Sequence Alignment Editor and Viewer (Version 2.8.2)
3 * Copyright (C) 2014 The Jalview Authors
5 * This file is part of Jalview.
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.
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.
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.
23 import jalview.datamodel.SequenceI;
26 * Assorted methods for analysing or comparing sequences.
28 public class Comparison
30 private static final int EIGHTY_FIVE = 85;
32 private static final int TO_UPPER_CASE = 'a' - 'A';
34 private static final char GAP_SPACE = ' ';
36 private static final char GAP_DOT = '.';
38 private static final char GAP_DASH = '-';
40 public static final String GapChars = new String(new char[]
41 { GAP_SPACE, GAP_DOT, GAP_DASH });
51 * @return DOCUMENT ME!
53 public static final float compare(SequenceI ii, SequenceI jj)
55 return Comparison.compare(ii, jj, 0, ii.getLength() - 1);
59 * this was supposed to be an ungapped pid calculation
71 public static float compare(SequenceI ii, SequenceI jj, int start, int end)
73 String si = ii.getSequenceAsString();
74 String sj = jj.getSequenceAsString();
76 int ilen = si.length() - 1;
77 int jlen = sj.length() - 1;
79 while (Comparison.isGap(si.charAt(start + ilen)))
84 while (Comparison.isGap(sj.charAt(start + jlen)))
95 for (int j = 0; j < jlen; j++)
97 if (si.substring(start + j, start + j + 1).equals(
98 sj.substring(start + j, start + j + 1)))
106 pid = (float) match / (float) ilen * 100;
110 for (int j = 0; j < jlen; j++)
112 if (si.substring(start + j, start + j + 1).equals(
113 sj.substring(start + j, start + j + 1)))
121 pid = (float) match / (float) jlen * 100;
128 * this is a gapped PID calculation
136 public final static float PID(String seq1, String seq2)
138 return PID(seq1, seq2, 0, seq1.length());
141 static final int caseShift = 'a' - 'A';
143 // Another pid with region specification
144 public final static float PID(String seq1, String seq2, int start, int end)
146 return PID(seq1, seq2, start, end, true, false);
150 * Calculate percent identity for a pair of sequences over a particular range,
151 * with different options for ignoring gaps.
160 * - if true - gaps match any character, if false, do not match
162 * @param ungappedOnly
163 * - if true - only count PID over ungapped columns
166 public final static float PID(String seq1, String seq2, int start,
167 int end, boolean wcGaps, boolean ungappedOnly)
169 int s1len = seq1.length();
170 int s2len = seq2.length();
172 int len = Math.min(s1len, s2len);
181 start = len - 1; // we just use a single residue for the difference
184 int elen = len - start, bad = 0;
188 for (int i = start; i < len; i++)
190 chr1 = seq1.charAt(i);
192 chr2 = seq2.charAt(i);
193 agap = isGap(chr1) || isGap(chr2);
194 if ('a' <= chr1 && chr1 <= 'z')
197 // Faster than toUpperCase
200 if ('a' <= chr2 && chr2 <= 'z')
203 // Faster than toUpperCase
231 return ((float) 100 * (elen - bad)) / elen;
235 * Answers true if the supplied character is a recognised gap character, else
236 * false. Currently hard-coded to recognise '-', '-' or ' ' (hyphen / dot /
243 public static final boolean isGap(char c)
245 return (c == GAP_DASH || c == GAP_DOT || c == GAP_SPACE) ? true : false;
249 * Answers true if more than 85% of the sequence residues (ignoring gaps) are
250 * A, G, C, T or U, else false. This is just a heuristic guess and may give a
251 * wrong answer (as AGCT are also animo acid codes).
256 public static final boolean isNucleotide(SequenceI[] seqs)
264 for (SequenceI seq : seqs)
266 for (char c : seq.getSequence())
268 if ('a' <= c && c <= 'z')
273 if (c == 'A' || c == 'G' || c == 'C' || c == 'T' || c == 'U')
277 else if (!Comparison.isGap(c))
285 * Check for nucleotide count > 85% of total count (in a form that evades
286 * int / float conversion or divide by zero).
288 if (ntCount * 100 > EIGHTY_FIVE * (ntCount + aaCount))