Context Navigation

source: branches/ralph/src/main/java/de/ugoe/cs/autoquest/tasktrees/alignment/algorithms/NeedlemanWunsch.java @ 1616

Last change on this file since 1616 was 1616, checked in by rkrimmel, 10 years ago
Small fixes
File size: 7.1 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats: