X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;ds=inline;f=src%2Fjalview%2Fanalysis%2FAlignSeq.java;h=832f5198781b804b055ee48f596c5bd6fc570b1c;hb=3889d8cc311503f8d40d9618acec1b0f4baf2aec;hp=c9ae0446605b4fa6ed1ac598379c789d6c76aed2;hpb=f55e45cc5d9db13675e9de77c45bd7b657da08ac;p=jalview.git diff --git a/src/jalview/analysis/AlignSeq.java b/src/jalview/analysis/AlignSeq.java index c9ae044..832f519 100755 --- a/src/jalview/analysis/AlignSeq.java +++ b/src/jalview/analysis/AlignSeq.java @@ -1,23 +1,21 @@ /* - * Jalview - A Sequence Alignment Editor and Viewer (Development Version 2.4.1) - * Copyright (C) 2009 AM Waterhouse, J Procter, G Barton, M Clamp, S Searle + * Jalview - A Sequence Alignment Editor and Viewer (Version 2.8) + * Copyright (C) 2012 J Procter, AM Waterhouse, LM Lui, J Engelhardt, 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 file is part of Jalview. * - * 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. + * 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 + * 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 + * You should have received a copy of the GNU General Public License along with Jalview. If not, see . */ package jalview.analysis; - import java.util.*; import java.awt.*; @@ -113,17 +111,17 @@ public class AlignSeq String type; - private int[] charToInt; + private int[] charToInt; /** * Creates a new AlignSeq object. * * @param s1 - * DOCUMENT ME! + * DOCUMENT ME! * @param s2 - * DOCUMENT ME! + * DOCUMENT ME! * @param type - * DOCUMENT ME! + * DOCUMENT ME! */ public AlignSeq(SequenceI s1, SequenceI s2, String type) { @@ -135,16 +133,16 @@ public class AlignSeq * Creates a new AlignSeq object. * * @param s1 - * DOCUMENT ME! + * DOCUMENT ME! * @param s2 - * DOCUMENT ME! + * DOCUMENT ME! * @param type - * DOCUMENT ME! + * DOCUMENT ME! */ public AlignSeq(SequenceI s1, String string1, SequenceI s2, String string2, String type) { - SeqInit(s1, string1, s2, string2, type); + SeqInit(s1, string1.toUpperCase(), s2, string2.toUpperCase(), type); } /** @@ -268,18 +266,13 @@ public class AlignSeq } /** - * DOCUMENT ME! - * - * @param s1 - * DOCUMENT ME! - * @param string1 - - * string to align for sequence1 - * @param s2 - * sequence 2 - * @param string2 - - * string to align for sequence2 + * Construct score matrix for sequences with standard DNA or PEPTIDE matrix + * @param s1 - sequence 1 + * @param string1 - string to use for s1 + * @param s2 - sequence 2 + * @param string2 - string to use for s2 * @param type - * DNA or PEPTIDE + * DNA or PEPTIDE */ public void SeqInit(SequenceI s1, String string1, SequenceI s2, String string2, String type) @@ -290,6 +283,14 @@ public class AlignSeq SeqInit(string1, string2); } + /** + * Construct score matrix for sequences with custom substitution matrix + * @param s1 - sequence 1 + * @param string1 - string to use for s1 + * @param s2 - sequence 2 + * @param string2 - string to use for s2 + * @param scoreMatrix - substitution matrix to use for alignment + */ public void SeqInit(SequenceI s1, String string1, SequenceI s2, String string2, ScoreMatrix scoreMatrix) { @@ -303,7 +304,7 @@ public class AlignSeq * construct score matrix for string1 and string2 (after removing any existing * gaps * - * @param string1 + * @param string1 * @param string2 */ private void SeqInit(String string1, String string2) @@ -383,13 +384,13 @@ public class AlignSeq { intToStr = pep; charToInt = ResidueProperties.aaIndex; - defInt = 23; + defInt = ResidueProperties.maxProteinIndex; } else if (type.equals(AlignSeq.DNA)) { intToStr = dna; charToInt = ResidueProperties.nucleotideIndex; - defInt = 4; + defInt = ResidueProperties.maxNucleotideIndex; } else { @@ -502,12 +503,25 @@ public class AlignSeq // TODO: Use original sequence characters rather than re-translated // characters in output // Find the biggest id length for formatting purposes + String s1id = s1.getName(), s2id = s2.getName(); int maxid = s1.getName().length(); if (s2.getName().length() > maxid) { maxid = s2.getName().length(); } - + if (maxid > 30) + { + maxid = 30; + // JAL-527 - truncate the sequence ids + if (s1.getName().length() > maxid) + { + s1id = s1.getName().substring(0, 30); + } + if (s2.getName().length() > maxid) + { + s2id = s2.getName().substring(0, 30); + } + } int len = 72 - maxid - 1; int nochunks = ((aseq1.length - count) / len) + 1; pid = 0; @@ -526,8 +540,7 @@ public class AlignSeq for (int j = 0; j < nochunks; j++) { // Print the first aligned sequence - output.append(new Format("%" + (maxid) + "s").form(s1.getName()) - + " "); + output.append(new Format("%" + (maxid) + "s").form(s1id) + " "); for (int i = 0; i < len; i++) { @@ -573,8 +586,7 @@ public class AlignSeq // Now print the second aligned sequence output = output.append("\n"); - output = output.append(new Format("%" + (maxid) + "s").form(s2 - .getName()) + output = output.append(new Format("%" + (maxid) + "s").form(s2id) + " "); for (int i = 0; i < len; i++) @@ -604,7 +616,7 @@ public class AlignSeq * DOCUMENT ME! * * @param mat - * DOCUMENT ME! + * DOCUMENT ME! */ public void printScoreMatrix(int[][] mat) { @@ -644,9 +656,9 @@ public class AlignSeq * DOCUMENT ME! * * @param i - * DOCUMENT ME! + * DOCUMENT ME! * @param j - * DOCUMENT ME! + * DOCUMENT ME! * * @return DOCUMENT ME! */ @@ -718,7 +730,7 @@ public class AlignSeq { E[i][0] = -gapOpen; F[i][0] = max(score[i - 1][0] - gapOpen, F[i - 1][0] - gapExtend); - + score[i][0] = max(lookup[seq1[i]][seq2[0]] * 10, E[i][0], F[i][0]); traceback[i][0] = -1; } @@ -742,9 +754,9 @@ public class AlignSeq * DOCUMENT ME! * * @param gapChar - * DOCUMENT ME! + * DOCUMENT ME! * @param seq - * DOCUMENT ME! + * DOCUMENT ME! * * @return DOCUMENT ME! */ @@ -765,11 +777,11 @@ public class AlignSeq * DOCUMENT ME! * * @param i1 - * DOCUMENT ME! + * DOCUMENT ME! * @param i2 - * DOCUMENT ME! + * DOCUMENT ME! * @param i3 - * DOCUMENT ME! + * DOCUMENT ME! * * @return DOCUMENT ME! */ @@ -794,9 +806,9 @@ public class AlignSeq * DOCUMENT ME! * * @param i1 - * DOCUMENT ME! + * DOCUMENT ME! * @param i2 - * DOCUMENT ME! + * DOCUMENT ME! * * @return DOCUMENT ME! */ @@ -816,9 +828,9 @@ public class AlignSeq * DOCUMENT ME! * * @param s - * DOCUMENT ME! + * DOCUMENT ME! * @param type - * DOCUMENT ME! + * DOCUMENT ME! * * @return DOCUMENT ME! */ @@ -839,11 +851,12 @@ public class AlignSeq try { seq1[i] = charToInt[c]; // set accordingly from setType - if (seq1[i]<0 || seq1[i] > defInt) // set from setType: 23 for peptides, or 4 for NA. + if (seq1[i] < 0 || seq1[i] > defInt) // set from setType: 23 for + // peptides, or 4 for NA. { seq1[i] = defInt; } - + } catch (Exception e) { seq1[i] = defInt; @@ -857,15 +870,15 @@ public class AlignSeq * DOCUMENT ME! * * @param g - * DOCUMENT ME! + * DOCUMENT ME! * @param mat - * DOCUMENT ME! + * DOCUMENT ME! * @param n - * DOCUMENT ME! + * DOCUMENT ME! * @param m - * DOCUMENT ME! + * DOCUMENT ME! * @param psize - * DOCUMENT ME! + * DOCUMENT ME! */ public static void displayMatrix(Graphics g, int[][] mat, int n, int m, int psize) @@ -907,4 +920,185 @@ public class AlignSeq } } } + + /** + * Compute a globally optimal needleman and wunsch alignment between two + * sequences + * + * @param s1 + * @param s2 + * @param type + * AlignSeq.DNA or AlignSeq.PEP + */ + public static AlignSeq doGlobalNWAlignment(SequenceI s1, SequenceI s2, + String type) + { + AlignSeq as = new AlignSeq(s1, s2, type); + + as.calcScoreMatrix(); + as.traceAlignment(); + return as; + } + + /** + * + * @return mapping from positions in S1 to corresponding positions in S2 + */ + public jalview.datamodel.Mapping getMappingFromS1(boolean allowmismatch) + { + ArrayList as1 = new ArrayList(), as2 = new ArrayList(); + int pdbpos = s2.getStart() + getSeq2Start() - 2; + int alignpos = s1.getStart() + getSeq1Start() - 2; + int lp2 = pdbpos - 3, lp1 = alignpos - 3; + boolean lastmatch = false; + // and now trace the alignment onto the atom set. + for (int i = 0; i < astr1.length(); i++) + { + char c1 = astr1.charAt(i), c2 = astr2.charAt(i); + if (c1 != '-') + { + alignpos++; + } + + if (c2 != '-') + { + pdbpos++; + } + + if (allowmismatch || c1 == c2) + { + lastmatch = true; + // extend mapping interval. + if (lp1 + 1 != alignpos) + { + as1.add(Integer.valueOf(alignpos)); + lp1 = alignpos; + } + + if (lp2 + 1 != pdbpos) + { + as2.add(Integer.valueOf(pdbpos)); + lp2 = pdbpos; + } + } + else + { + lastmatch = false; + } + } + // construct range pairs + int[] mapseq1 = new int[as1.size() + (lastmatch ? 1 : 0)], mapseq2 = new int[as2 + .size() + (lastmatch ? 1 : 0)]; + int i = 0; + for (Integer ip : as1) + { + mapseq1[i++] = ip; + } + ; + i = 0; + for (Integer ip : as2) + { + mapseq2[i++] = ip; + } + ; + if (lastmatch) + { + mapseq1[mapseq1.length - 1] = alignpos; + mapseq2[mapseq2.length - 1] = pdbpos; + } + MapList map = new MapList(mapseq1, mapseq2, 1, 1); + + jalview.datamodel.Mapping mapping = new Mapping(map); + mapping.setTo(s2); + return mapping; + } + + /** + * compute the PID vector used by the redundancy filter. + * + * @param originalSequences + * - sequences in alignment that are to filtered + * @param omitHidden + * - null or strings to be analysed (typically, visible portion of + * each sequence in alignment) + * @param start + * - first column in window for calculation + * @param end + * - last column in window for calculation + * @param ungapped + * - if true then use ungapped sequence to compute PID + * @return vector containing maximum PID for i-th sequence and any sequences + * longer than that seuqence + */ + public static float[] computeRedundancyMatrix( + SequenceI[] originalSequences, String[] omitHidden, int start, + int end, boolean ungapped) + { + int height = originalSequences.length; + float[] redundancy = new float[height]; + int[] lngth = new int[height]; + for (int i = 0; i < height; i++) + { + redundancy[i] = 0f; + lngth[i] = -1; + } + + // long start = System.currentTimeMillis(); + + float pid; + String seqi, seqj; + for (int i = 0; i < height; i++) + { + + for (int j = 0; j < i; j++) + { + if (i == j) + { + continue; + } + + if (omitHidden == null) + { + seqi = originalSequences[i].getSequenceAsString(start, end); + seqj = originalSequences[j].getSequenceAsString(start, end); + } + else + { + seqi = omitHidden[i]; + seqj = omitHidden[j]; + } + if (lngth[i] == -1) + { + String ug = AlignSeq.extractGaps(Comparison.GapChars, seqi); + lngth[i] = ug.length(); + if (ungapped) + { + seqi = ug; + } + } + if (lngth[j] == -1) + { + String ug = AlignSeq.extractGaps(Comparison.GapChars, seqj); + lngth[j] = ug.length(); + if (ungapped) + { + seqj = ug; + } + } + pid = Comparison.PID(seqi, seqj); + + // use real sequence length rather than string length + if (lngth[j] < lngth[i]) + { + redundancy[j] = Math.max(pid, redundancy[j]); + } + else + { + redundancy[i] = Math.max(pid, redundancy[i]); + } + + } + } + return redundancy; + } }