source: branches/autoquest-core-tasktrees-alignment/src/main/java/de/ugoe/cs/autoquest/tasktrees/alignment/algorithms/SmithWaterman.java @ 1734

Last change on this file since 1734 was 1734, checked in by rkrimmel, 10 years ago

Added automatically created javadoc, still needs to be commented properly though
File size: 8.1 KB
Line 
1/*
2 *
3 */
4package de.ugoe.cs.autoquest.tasktrees.alignment.algorithms;
5
6import java.util.ArrayList;
7import java.util.logging.Level;
8
9import de.ugoe.cs.autoquest.tasktrees.alignment.matrix.SubstitutionMatrix;
10import de.ugoe.cs.autoquest.tasktrees.alignment.algorithms.Constants;
11import de.ugoe.cs.util.console.Console;
12
13// TODO: Auto-generated Javadoc
14/**
15 * The Class SmithWaterman.
16 */
17public class SmithWaterman implements AlignmentAlgorithm {
18
19        /** The first input. */
20        private int[] input1;
21
22        /** The second input String. */
23        private int[] input2;
24
25        /** The lengths of the input. */
26        private int length1, length2;
27
28        /**
29         * The score matrix. The true scores should be divided by the normalization
30         * factor.
31         */
32        private MatrixEntry[][] matrix;
33
34        /** The alignment. */
35        private ArrayList<NumberSequence> alignment;
36
37        /** Substitution matrix to calculate scores. */
38        private SubstitutionMatrix submat;
39
40        /* (non-Javadoc)
41         * @see de.ugoe.cs.autoquest.tasktrees.alignment.algorithms.AlignmentAlgorithm#align(de.ugoe.cs.autoquest.tasktrees.alignment.algorithms.NumberSequence, de.ugoe.cs.autoquest.tasktrees.alignment.algorithms.NumberSequence, de.ugoe.cs.autoquest.tasktrees.alignment.matrix.SubstitutionMatrix, float)
42         */
43        @Override
44        public void align(NumberSequence input1, NumberSequence input2,
45                        SubstitutionMatrix submat, float threshold) {
46                this.input1 = input1.getSequence();
47                this.input2 = input2.getSequence();
48                length1 = input1.size();
49                length2 = input2.size();
50                this.submat = submat;
51
52                // System.out.println("Starting SmithWaterman algorithm with a "
53                // + submat.getClass() + " Substitution Matrix: " +
54                // submat.getClass().getCanonicalName());
55
56                matrix = new MatrixEntry[length1 + 1][length2 + 1];
57                alignment = new ArrayList<NumberSequence>();
58
59                for (int i = 0; i < (length1 + 1); i++) {
60                        for (int j = 0; j < (length2 + 1); j++) {
61                                matrix[i][j] = new MatrixEntry();
62                        }
63                }
64
65                buildMatrix();
66                traceback();
67
68        }
69
70        /**
71         * Build the score matrix using dynamic programming.
72         */
73        private void buildMatrix() {
74                if (submat.getGapPenalty() >= 0) {
75                        throw new Error("Indel score must be negative");
76                }
77
78                // it's a gap
79                matrix[0][0].setScore(0);
80                matrix[0][0].setPrevious(null); // starting point
81
82                // the first column
83                for (int j = 1; j < length2; j++) {
84                        matrix[0][j].setScore(0);
85                        matrix[0][j].setPrevious(matrix[0][j - 1]);
86                        matrix[0][j].setYvalue(input2[j]);
87                        matrix[0][j].setXvalue(Constants.GAP_SYMBOL);
88                }
89                // the first row
90
91                for (int j = 1; j < length1; j++) {
92                        matrix[j][0].setScore(0);
93                        matrix[j][0].setPrevious(matrix[j - 1][0]);
94                        matrix[j][0].setXvalue(input1[j]);
95                        matrix[j][0].setYvalue(Constants.GAP_SYMBOL);
96                }
97
98                for (int i = 1; i < length1; i++) {
99
100                        for (int j = 1; j < length2; j++) {
101                                final double diagScore = matrix[i - 1][j - 1].getScore()
102                                                + similarity(i, j);
103                                final double upScore = matrix[i][j - 1].getScore()
104                                                + submat.getGapPenalty();
105                                final double leftScore = matrix[i - 1][j].getScore()
106                                                + submat.getGapPenalty();
107
108                                matrix[i][j].setScore(Math.max(diagScore,
109                                                Math.max(upScore, Math.max(leftScore, 0))));
110
111                                // find the directions that give the maximum scores.
112                                // TODO: Multiple directions are ignored, we choose the first
113                                // maximum score
114                                // True if we had a match
115                                if (diagScore == matrix[i][j].getScore()) {
116                                        matrix[i][j].setPrevious(matrix[i - 1][j - 1]);
117                                        matrix[i][j].setXvalue(input1[i - 1]);
118                                        matrix[i][j].setYvalue(input2[j - 1]);
119                                }
120                                // true if we took an event from sequence x and not from y
121                                if (leftScore == matrix[i][j].getScore()) {
122                                        matrix[i][j].setXvalue(input1[i - 1]);
123                                        matrix[i][j].setYvalue(Constants.GAP_SYMBOL);
124                                        matrix[i][j].setPrevious(matrix[i - 1][j]);
125                                }
126                                // true if we took an event from sequence y and not from x
127                                if (upScore == matrix[i][j].getScore()) {
128                                        matrix[i][j].setXvalue(Constants.GAP_SYMBOL);
129                                        matrix[i][j].setYvalue(input2[j - 1]);
130                                        matrix[i][j].setPrevious(matrix[i][j - 1]);
131                                }
132                                if (0 == matrix[i][j].getScore()) {
133                                        matrix[i][j].setPrevious(matrix[0][0]);
134                                        matrix[i][j].setXvalue(-2);
135                                        matrix[i][j].setYvalue(-2);
136                                }
137                        }
138
139                        // Set the complete score cell
140
141                }
142        }
143
144        /*
145         * (non-Javadoc)
146         *
147         * @see
148         * de.ugoe.cs.autoquest.tasktrees.alignment.algorithms.AlignmentAlgorithm
149         * #getAlignment()
150         */
151        @Override
152        public ArrayList<NumberSequence> getAlignment() {
153                return alignment;
154        }
155
156        /*
157         * (non-Javadoc)
158         *
159         * @see
160         * de.ugoe.cs.autoquest.tasktrees.alignment.algorithms.AlignmentAlgorithm
161         * #getAlignmentScore()
162         */
163        @Override
164        public double getAlignmentScore() {
165                return getMaxScore();
166        }
167
168        /* (non-Javadoc)
169         * @see de.ugoe.cs.autoquest.tasktrees.alignment.algorithms.AlignmentAlgorithm#getMatches()
170         */
171        @Override
172        public ArrayList<Match> getMatches() {
173                // TODO Auto-generated method stub
174                return null;
175        }
176
177        /**
178         * Get the maximum value in the score matrix.
179         *
180         * @return the max score
181         */
182        @Override
183        public double getMaxScore() {
184                double maxScore = 0;
185
186                // skip the first row and column
187                for (int i = 1; i <= length1; i++) {
188                        for (int j = 1; j < length2; j++) {
189                                if (matrix[i][j].getScore() > maxScore) {
190                                        maxScore = matrix[i][j].getScore();
191                                }
192                        }
193                }
194
195                return maxScore;
196        }
197
198        /* (non-Javadoc)
199         * @see de.ugoe.cs.autoquest.tasktrees.alignment.algorithms.AlignmentAlgorithm#printAlignment()
200         */
201        @Override
202        public void printAlignment() {
203                MatrixEntry tmp = matrix[length1 + 1][0];
204                String aligned1 = "";
205                String aligned2 = "";
206                int count = 0;
207                do {
208                        String append1 = "";
209                        String append2 = "";
210
211                        if (tmp.getXvalue() == Constants.GAP_SYMBOL) {
212                                append1 = "  ___";
213                        } else {
214                                append1 = String.format("%5d", tmp.getXvalue());
215                        }
216
217                        if (tmp.getYvalue() == Constants.GAP_SYMBOL) {
218                                append2 = "  ___";
219                        } else {
220                                append2 = String.format("%5d", tmp.getYvalue());
221                        }
222                        if (count != 0) {
223                                aligned1 = append1 + aligned1;
224                                aligned2 = append2 + aligned2;
225                        }
226
227                        tmp = tmp.getPrevious();
228                        count++;
229
230                } while (tmp != null);
231                System.out.println(aligned1);
232                System.out.println(aligned2);
233        }
234
235        /**
236         * print the dynmaic programming matrix.
237         */
238        @Override
239        public void printDPMatrix() {
240                System.out.print("          ");
241                for (int i = 1; i <= length1; i++) {
242                        System.out.format("%5d", input1[i - 1]);
243                }
244                System.out.println();
245                for (int j = 0; j <= length2; j++) {
246                        if (j > 0) {
247                                System.out.format("%5d ", input2[j - 1]);
248                        } else {
249                                System.out.print("      ");
250                        }
251                        for (int i = 0; i <= length1; i++) {
252                                System.out.format("%4.1f ", matrix[i][j].getScore());
253                        }
254                        System.out.println();
255                }
256        }
257
258        /**
259         * Sets the alignment.
260         *
261         * @param alignment the new alignment
262         */
263        public void setAlignment(ArrayList<NumberSequence> alignment) {
264                this.alignment = alignment;
265        }
266
267        /**
268         * Compute the similarity score of substitution The position of the first
269         * character is 1. A position of 0 represents a gap.
270         *
271         * @param i
272         *            Position of the character in str1
273         * @param j
274         *            Position of the character in str2
275         * @return Cost of substitution of the character in str1 by the one in str2
276         */
277        private double similarity(int i, int j) {
278                return submat.getScore(input1[i - 1], input2[j - 1]);
279        }
280
281        /**
282         * Traceback.
283         */
284        public void traceback() {
285                MatrixEntry tmp = matrix[length1][length2];
286                final int aligned1[] = new int[length1 + length2 + 2];
287                final int aligned2[] = new int[length1 + length2 + 2];
288                int count = 0;
289                do {
290                        if ((length1 + length2 + 2) == count) {
291                                Console.traceln(Level.WARNING,
292                                                "Traceback longer than both sequences summed up!");
293                                break;
294                        }
295                        aligned1[count] = tmp.getXvalue();
296                        aligned2[count] = tmp.getYvalue();
297
298                        tmp = tmp.getPrevious();
299                        count++;
300
301                } while (tmp != null);
302                count--;
303                // reverse order of the alignment
304                final int reversed1[] = new int[count];
305                final int reversed2[] = new int[count];
306
307                for (int i = count - 1; i > 0; i--) {
308                        reversed1[reversed1.length - i] = aligned1[i];
309                        reversed2[reversed2.length - i] = aligned2[i];
310                }
311
312                final NumberSequence ns1 = new NumberSequence(reversed1.length);
313                final NumberSequence ns2 = new NumberSequence(reversed2.length);
314                ns1.setSequence(reversed1);
315                ns2.setSequence(reversed2);
316
317                alignment.add(ns1);
318                alignment.add(ns2);
319        }
320
321}
Note: See TracBrowser for help on using the repository browser.