Wstęp

  • Programowanie dynamiczne to algorytmiczny sposób rozwiązywania problemów
  • Opiera się na podziale problemu na mniejsze części łatwiejsze do rozwiązania
  • W przeciwieństwie do metody dziel-i-rządź mniejsze podproblemy nie są rozłączne

Algorytm Needlemana-Wunscha

  • Służy do znalezienia uliniowienia (ang. alignment) różnej długości sekwencji s i t
  • Rezultatem działania algorytmu są sekwencje s' oraz t' , które spełniają następujące cechy:
    • Są równej długości
    • Oprócz alfabetu używanego w oryginalnych sekwencjach s i t mogą zawierać specjalny znak pusty (ang. gap) oznaczany myślnikiem '-'
    • Na każdej pozycji s' i t' para znaków może być:
      • równa np. AA (ang. match )
      • różna np. AU (ang. mismatch )
      • zawierająca znak pusty np. -U (ang. indel , insertion/deletion)
  • Istnieją różne schematy punktowania poszczególnych sytuacji; dla uproszczenia przyjmijmy, że match = +1, mismatch = -1, indel = -1
  • Programowanie dynamiczne polega na tym, żeby dla coraz większych podciągów s i t podejmować decyzje o przyznawaniu match , mismatch lub indel

Algorytm krok po kroku:

  • Przygotuj macierz w taki sposób, żeby każdy znak sekwencji s był w osobnej kolumnie, a każdy znak t w osobnym wierszu; dodaj na początku jeden dodatkowy wierz i kolumnę
  PRETTY
        
P       
R       
T       
T       
E       
I       
N       
  • Wypełnij dodatkowy wiersz i kolumnę malejącymi wartościami
  PRETTY
 0-1-2-3-4-5-6
P-1      
R-2      
T-3      
T-4      
E-5      
I-6      
N-7      
  • Od teraz macierz można wypełniać tylko dla tych pól, które mają znaną wartość dla sąsiada z lewej, z góry i po ukosie między nimi
  • Każde pole macierzy wypełniamy wartością maksymalną z:
    • Komórki po lewej minus 1 ( indel )
    • Komórki na górze minus 1 ( indel )
    • Komórki po ukosie minus 1 jeżeli jest mismatch
    • Komórki po ukosie plus 1 jeżeli jest match
  • Zatem wypełniamy cały pierwszy wiersz. Na przecięciu P i P pojawia się 1, bo mamy match , a po ukosie wartość 0. Można zauważyć, że wszystkie dalsze wartości w wierszu stosują krok "komórka po lewej minus 1 (indel)"
  PRETTY
 0-1-2-3-4-5-6
P-110-1-2-3-4
R-2      
T-3      
T-4      
E-5      
I-6      
N-7      
  • Analogicznie wypełniamy pozostałą część macierzy
  PRETTY
 0-1-2-3-4-5-6
P-110-1-2-3-4
R-20210-1-2
T-3-111210
T-4-200232
E-5-3-11122
I-6-4-20011
N-7-5-3-1-100
  • Kluczowe jest teraz odczytanie wyniku z macierzy. Musimy prześledzić wstecz (ang. traceback) skąd wzięła się każda wartość zaczynając od ostatniego pola w macierzy w prawym dolnym rogu
  • Śledząc wstecz musimy sprawdzić czy aktualna wartość powstała z operacji indel od komórki z góry lub operacji indel od komórki z lewej, a może przez match lub mismatch
  • Czasami jest to niejednoznaczne. Już w pierwszym kroku wstecz nie wiadomo czy wybrać  indel z komórki u góry czy mismatch z komórki po ukosie. W praktyce można przyjąć stałą kolejność sprawdzanych wyborów np. najpierw indel z góry, potem indel z lewej i na końcu mismatch
  • Macierz po wykonaniu śledzenia wstecz wygląda następująco
  PRETTY
 0-1-2-3-4-5-6
P-110-1-2-3-4
R-20210-1-2
T-3-111210
T-4-200232
E-5-3-11122
I-6-4-20011
N-7-5-3-1-100
  • Aby wypisać wynik należy przy każdym kroku śledzenia tworzyć sekwencje s' oraz t' :
    1. W górę, zatem s' = -, t' = N
    2. W górę, zatem s' = --, t' = IN
    3. Na ukos, zatem s' = Y--, t' = EIN
    4. Na ukos, zatem s' = TY--, t' = TEIN
    5. Na ukos, zatem s' = TTY--, t' = TTEIN
    6. W lewo, zatem s' = ETTY--, t' = -TTEIN
    7. Na ukos, zatem s' = RETTY--, t' = R-TTEIN
    8. Na ukos, zatem s' = PRETTY--, t' = PR-TTEIN
  • Ostateczny alignment
PRETTY--
PR-TTEIN

Napisz program, który na standardowym wejściu otrzyma dwie linie (każda z sekwencją) a na standardowe wyjście wypisze alignment wyznaczony algorytmem Needlemana-Wunscha.

Przykład:

$ echo -e 'PRETTY\nPRTTEIN' | ./needleman-wunsch.py
PRETTY--
PR-TTEIN

Algorytm Nussinov

  • W 1980 Ruth Nussinov podała algorytm przewidujący strukturę drugorzędową RNA
  • Algorytm ten wykorzystuje programowanie dynamiczne by dla zadanej sekwencji znaleźć maksymalną liczbę par kanonicznych, które da się przypisać

Algorytm krok po kroku

  • Budujemy macierz kwadratową gdzie wiersze i kolumny reprezentują pojedyncze znaki sekwencji
 GAUUACA
G       
A       
U       
U       
A       
C       
A       
  • Wypełniamy przekątną oraz jej przesunięcie w dół zerami
 GAUUACA
G0      
A00     
U 00    
U  00   
A   00  
C    00 
A     00
  • W tym algorytmie wypełniamy macierz podążając od przekątnej aż do prawego górnego rogu
  • Dla każdej komórki M[i,j] macierzy wyznaczamy maksimum z:
    • M[i + 1][j - 1] + s(i, j) , gdzie s(i, j) to 1 gdy na przecięciu macierzy mamy parę kanoniczną G-C lub A-U i 0 w innym wypadku
    • M[i][k] + M[k + 1][j] dla każdego i ≤ k < j
  • Wizualnie prezentuje się to tak, że aby wyznaczyć ostateczny wynik, musimy albo wybrać z albo sumy pól oznaczonych jako t lub u lub v lub w lub x lub y. Oczywiście w toku rozwiązywania problemu wszystkie puste pola macierzy zostaną wypełnione.
 GAUUACA
Gt=0uvwxy?
A00   zt
U 00   u
U  00  v
A   00 w
C    00x
A     0y
  • Wypełniamy zatem pierwszą serię pól. Na początku tylko podstawowe dwie pary A-U możemy dopisać
 GAUUACA
G00     
A001    
U 000   
U  001  
A   000 
C    000
A     00
  • Kolejna iteracja. W pierwszym wierszu ( i = 1, j = 3 ) mamy G-U, którego nie zaliczamy jako pary, ale dla k = 1 , mamy M[1,1] + M[2,3] = 0 + 1 = 1 , dlatego w macierzy wpisujemy 1
 GAUUACA
G001    
A0011   
U 0001  
U  0011 
A   0000
C    000
A     00
  • Kolejna iteracja. W drugim wierszu ( i = 2, j = 5 ) mamy wartość 2, ponieważ zachodzi M[2][3] + M[4][5] = 1 + 1 = 2
 GAUUACA
G0011   
A00112  
U 00011 
U  00111
A   0000
C    000
A     00
  • Kolejna iteracja. W pierwszym wierszu ( i = 1, j = 5 ) pojawia się 2, ponieważ dla k = 3 zachodzi M[1][3] + M[4][5] = 1 + 1 = 2
 GAUUACA
G00112  
A001122 
U 000112
U  00111
A   0000
C    000
A     00
  • Kolejna iteracja
 GAUUACA
G001123 
A0011222
U 000112
U  00111
A   0000
C    000
A     00
  • Ostatnia iteracja
 GAUUACA
G0011233
A0011222
U 000112
U  00111
A   0000
C    000
A     00
  • Wynik odczytujemy z lewego górnego narożnika. Jest to maksymalna liczba par, jakie mogą się znaleźć w strukturze o zadanej sekwencji
  • Śledzenie wsteczne ponownie wymaga by na każdym kroku odczytać skąd wzięła się dana wartość w macierzy
    1. Wartość  M[1][7] to M[1][6] + M[7][7] (3 + 0)
    2. Wartość  M[1][6] to M[2][5] + 1 (para G-C)
    3. Wartość M[2][5] to M[2][3] + M[4][5] (1 + 1)
    4. Wartość M[2][3] to M[3][2] + 1 (para A-U)
    5. Wartość  M[4][5] to M[5][4] + 1 (para A-U)
  • A zatem zestaw par to (1,6), (2,3) i (4,5), a struktura dot-bracket:
     
    GAUUACA
    (()()). 

Napisz program, który odczyta ze standardowego wejścia sekwencję, a wypisze strukturę w postaci dot-bracket np.

$ echo GAUUACA | ./nussinov.py 
GAUUACA
(()()). 

  • W obu zadaniach jeżeli jest niejednoznaczność to można wybrać dowolnie. Będę sprawdzał rozwiązania pod kątem tego czy są optymalne. Np. dla UAA prawidłowe będzie zarówno (). oraz (.)
  • Zadania są traktowane niezależnie tym razem, więc za każde jest do uzyskania 5 punktów