1 package fr.orsay.lri.varna.applications;
3 import java.util.ArrayList;
4 import java.util.Hashtable;
6 import fr.orsay.lri.varna.models.rna.ModeleBP;
8 public class SecStrConsensus {
12 * Internal class to represent a simple base-pair.
16 static class SimpleBP{
20 public SimpleBP(int i5, int i3)
27 public static int[] extractConsensus(ArrayList<ArrayList<SimpleBP>> bps)
29 // We do not currently know the length of the alignment
30 // => Estimate it as the biggest index of a base-pair plus one.
32 for (ArrayList<SimpleBP> strs : bps)
34 for (SimpleBP bp : strs)
36 maxlength = Math.max(1+Math.max(bp.bp5, bp.bp3), maxlength);
39 // Now we have a good estimate for length, allocate and initialize data
40 // to be fed to the dynamic programming procedure.
41 ArrayList<Hashtable<Integer,Double>> seq = new ArrayList<Hashtable<Integer,Double>>();
42 for (int i=0;i<seq.size();i++)
43 { seq.add(new Hashtable<Integer,Double>()); }
44 for (ArrayList<SimpleBP> strs : bps)
46 for (SimpleBP bp : strs)
50 Hashtable<Integer,Double> h = seq.get(i);
51 if (!h.containsKey(j))
55 h.put(j, h.get(i)+1.);
58 // At this point, seq contains, at each position i, a hashtable which associates,
59 // to each possible end j, the number of time a base-pair (i,j) occurs in the alignment
61 // We can now run the dynamic programming procedure on this data
62 double[][] mat = fillMatrix(seq);
63 ArrayList<SimpleBP> res = backtrack(mat,seq);
65 // Convert it to an array, ie finalres[i] = j >= 0 iff a base-pair (i,j) is present
66 // in the consensus, or -1 otherwise
67 int[] finalres = new int[seq.size()];
68 for (int i=0;i<seq.size();i++)
70 for (SimpleBP bp : res)
72 finalres[bp.bp5] = bp.bp3;
73 finalres[bp.bp3] = bp.bp5;
79 private static boolean canBasePair(ArrayList<Hashtable<Integer,Double>> seq, int i, int k)
81 return seq.get(i).containsKey(k);
84 // Returns the score of a potential base-pair, ie the number of structures in which it is found.
85 private static double basePairScore(ArrayList<Hashtable<Integer,Double>> seq, int i, int k)
87 return seq.get(i).get(k);
91 private static double[][] fillMatrix(ArrayList<Hashtable<Integer,Double>> seq)
94 double[][] tab = new double[n][n];
97 for(int i=0;i<n-m+1;i++)
103 tab[i][j] = Math.max(tab[i][j], tab[i+1][j]);
104 for (int k=i+1;k<=j;k++)
106 if (canBasePair(seq,i,k))
111 fact1 = tab[i+1][k-1];
118 tab[i][j] = Math.max(tab[i][j],basePairScore(seq,i,k)+fact1+fact2);
127 private static ArrayList<SimpleBP> backtrack(double[][] tab,ArrayList<Hashtable<Integer,Double>> seq)
129 return backtrack(tab,seq,0,seq.size()-1);
132 private static ArrayList<SimpleBP> backtrack(double[][] tab,ArrayList<Hashtable<Integer,Double>> seq, int i, int j)
134 ArrayList<SimpleBP> result = new ArrayList<SimpleBP>();
137 ArrayList<Integer> indices = new ArrayList<Integer>();
139 for (int k=i+1;k<=j;k++)
143 for (int k : indices)
147 if (tab[i][j] == tab[i+1][j])
149 result = backtrack(tab, seq, i+1,j);
154 if (canBasePair(seq,i,k))
159 fact1 = tab[i+1][k-1];
166 if (tab[i][j]==basePairScore(seq,i,k)+fact1+fact2)
168 result = backtrack(tab, seq, i+1,k-1);
169 result.addAll(backtrack(tab, seq, k+1,j));
170 result.add(new SimpleBP(i,k));