Category: WBO

Dziś rozmawialiśmy o algorytmach przybliżonego wyszukiwania sekwencji w bazach danych. Slajdy są tu: wyk8

Na ćwiczeniach spróbujemy zająć się wykorzystaniem programu blast:

0. Wczytaj plik w formacie FastQ microbial_reads.fastq (jest już na naszym serwerze jupyter, nie trzeba go tam ładować), przy pomocy SeqIO.parse(). Są to odczyty z mikrobiomu jelitowego myszy.

1. Wybierz kilka losowych, dość długich sekwencji DNA i uruchom dla nich program BLAST online, obejrzyj wyniki (jeśli nic się nie “trafiło”, możesz wybrać inne, dłuższe sekwencje)

2. Wykonaj wyszukiwanie dla tych samych sekwencji przy pomocy interfejsu API biopython’a do blasta online (NCBIWWW) i parsera xml (NCBIXML)

3. Znajdź, która z twoich wybranych sekwencji ma najistotniejsze trafienia do bazy NCBI

4. Przeanalizuj “trafioną” sekwencję. Czy to możliwe, aby ta sekwencja była dokładnie z tego organizmu, który był badany? Porównaj z wynikami algorytmu Smith’a-Waterman’a

 

 

Dzisiaj mówiliśmy o zastosowaniach HMM do reprezentacji profilów uliniowień oraz wyszukiwania genów.

Slajdy są tu: wyk7

W ramach ćwiczeń, warto jest zapoznać się z bazą PFAM  i możliwością wyszukiwania profilów domen przy pomocy programu HMMER. Mieliśmy też wykorzystać moduł HmmerIO, ale niestety nie działa on po niedawnych zmianach formatu HMMER.

Zasadniczo powinniśmy móc wykorzystać to narzędzie (pfam+hmmer) do wyszukania domen. I znalezienia ich modeli. Proponuję zacząć od białka tramtrack , które posiada 2 istotnie różne wersje TTK69 i TTK88. Proszę znaleźć te dwie sekwencje aminokwasowe w bazie flybase, ściągnąć je do siebie i wrzucić do wyszukiwarki domen pfam/hmmer. Następnie porównać struktury domen tych białek i zobaczyć jak to się ma do struktury egzonów genu tramtrack. Warto zobaczyć jak wyglądają opisy domen w postaci profilów hmm.

Dziś na wykładzie omówiliśmy ukryte modele Markowa i algorytmy rekonstrukcji parametrów z danych. Slajdy są tu: wyk6

Przy okazji – na naszym serwerze jupyter doinstalowalem mozliwosc uzywania python’a 3, ale niestety przestala dzialac funkcjonalnosc python’a 2. Jesli ktos chce korzystac ze stworzonych wczesniej notebooków, to trzeba wybrac opcję: kernel->change kernel-> python 3

Zadania na dziś:

0. Zapoznaj się z modułem hmmlearn. Będziemy go wykorzystywać do uczenia modeli z emisjami

1. Wyspy CpG znajdują się często w genomach, w szczególności genomie ludzkim. Spróbuj zdefiniować ukryty model Markowa, który ma 2 stany i spróbuj nauczyć go na sekwncji zawartej w pliku cpg.fa. Zrób to zarówno dla sekwencji liter (ACGT), jak i dla sekwencji dwunukleotydów (AA,AC,AG,AT, itp…) Czy możesz zinterpretować macierze emisji i przypisać jeden ze stanów do wysp CpG? Wykonaj kilkakrotnie proces uczenia (Baum-Welch) i zobacz czy wyniki są podobne. Jak interpretujesz prawdopodobieństwa w macierzy przejść. Czy coś możesz powiedzieć o średniej długości wysp CpG?

Warto przyjrzeć się przykładowi użycia klasy MultinomialHMM

2. czasami potrzebujemy użyć łańcuchów Markowa do segmentacji sygnału. Weżmy przykładowy plik – dane o methylacji histonów w rzodkiewniku (H3K9me2-crh6-3-chr1.sgr) zawierający poziom metylacji w różnych pozycjach (plik jest tekstowy, każda linia zawiera informacje o pozycji i wartośći). Wykorzystaj HMM z emisjami Gaussowskimi do segmentacji genomu na 2 lub 3 stany. Jakie są wysetymowane wartości średnie dla różnych stanów? Przykład wykorzystania HMMów z emisjami Gaussowskimi można znaleźć tu 

praca domowa:

Wykorzystaj model nauczony na danych o CpG i przetestuj które z 30 sekwencji w pliku cpg_test.fa są naprawdę wyspami CpG. Jako wynik proszę przysłać program w pythonie i wynik w pliku tekstowym.

Naszym zadaniem jest elegancka implementacja progresywnego uliniowienia wg następujących wytycznych

1) Zaimplementuj klasę profile, która zachowuje się podobnie do multiple alignment, ale kiedy żądamy indeksu a[i], to dostajemy i-tą kolumnę uliniowienia, zamiast wiersza (w Multiple Alignment z Biopythona mamy a[i] = wiersz, a[:,i]=kolumna)  (2 pkt)

2) zaimplementuj funkcje callback dla metody Bio.pairwise2.globalcc dla afinicznej kary za przerwę (parametryzowalnej 2 parametrami) i oceny substytucji profili, dla funkcji oceny sumy par (3 pkt)

3) zaimplementuj metodę obliczającą odległości sekwencji białkowych według prostego modelu estymacji odległości, gdzie d= -19/20 * ln (1 – (20/19)*p), dla p=procent zmienionych aminokwasów w globalnym uliniowieniu pary, oraz pozwalającą wyliczać macierz odległości parami (2 pkt)

4) Zaimplementuj metodę konstruującą uliniowienie progresywne korzystającą z drzewa stworzonego metodą neighbor joining (nj w Biopythonie) dla macierzy odległości uzyskanej metodą z p.3) i uliniowienie profili metodą z punktu 1) i 2). (3 pkt)

Wszystko to powinno być włożone do udokumentowanego modułu python’a i wysłane do mnie w emailu z “[WBO-1-2018]” w temacie do 16. IV 2018.

Dzisiaj mówiliśmy o metodach uliniowienia wielu sekwencji. Slajdy są tu: wyk5

Na zajęciach naszym celem jest zapoznanie się z metodami uliniowienia wielu sekwencji w pakiecie Biopython oraz metodami uliniowienia progresywnego. Warto przeczytać odpowiedni rozdział tutorialu biopythona

0. Tutaj warto, aby każdy z Państwa założył sobie na serwerze Jupyter własny folder, żeby pliki nam się nie mieszały. Proponuję numer indeksu, albo inny identyfikator jednoznaczny.

1.Weźmy na początek naszą rodzinę sekwencji białek histonowych. Jest dostępna na naszym serwerze jupyter w pliku histones.fa_ Na początek potrzebujemy wczytać je z pliku (SeqIO.parse ) przetłumaczyć je na sekwencje białkowe (metoda translate()) i zapisać do pliku fasta (SeqIO.write) wraz z odpowiednimi identyfikatorami i usuniętymi kodonami stopu (*, która powstaje w wyniku translate)

2. Wykorzystajmy program clustalw (Bio.Align.Applications.ClustalwCommandline) do wykonania mulituliniowienia tych sekwencji na naszym serwerze plik wykonywalny jest zainstalowany w katalogu /usr/bin. Po wykonaniu zadania wczytaj uliniowienie z pliku .aln (AlignIO.read format “clustal”) oraz drzewo filogenetyczne z pliku .dnd (Phylo.read format=”newick”)

3. Wykorzystaj teraz program muscle (Bio.Align.Applications.MuscleCommandline) do wykonania multiuliniowienia tych samych sekwencji. Wczytaj je z pliku w formacie fasta

4. Stwórz drzewo filogenetyczne na podstawie uzyskanego uliniowienia metodą neighbor joining jak na poprzednich zajęciach. Porównaj to drzewo do drzewa otrzymanego z programu clustalw

5. Powtórz ćwiczenia 2-4 dla większego zbioru sekwencji (np. Human_PAH_homologs). Czy tutaj widoczne są większe różnice między metodami?

6(*). Napisz prostą implementację uliniowienia progresywnego, która wykorzystuje globalne uliniowienie z modułu Bio.pairwise2 i funkcje “callback”

 

Dziś rozmawiamy o prostych metodach kostrukcji drzew na podstawie macierzy odległości. Slajdy na dzis wyk4

Na laboratorium będziemy konstruować drzewa w praktyce.

Warto zapoznać się z dokumentacją modułu Bio.Phylo w tutorialu jak i na stronie wiki.

1. Na początek ustalmy dwie listy sekwencji. Rozważmy sekwencje paralogiczne histonów z drożdży z poprzednich zajęć oraz ortologi ludzkiej hydroksylazy fenyloalaniny   Human_PAH_orthologues.fa. Dla histonów musimy je najpierw przetłumaczyć na białka (s.translate()), a dla drugiego zbioru mamy już sekwencje aminokwasowe.

2. Wylicz macierze odległości dla tych grup sekwencji przy pomocy klasy Bio.Phylo.TreeConstruction.DistanceCalculator  (dla macierzy BLOSUM62, osobno dla histonów, i osobno dla genów PAH – to  zajmie chwilę). Można użyć plików  Human_PAH_orthologues i Human_PAH_paralogues oraz Human_H2BFS_paralogues i duży plik Human_PAH_orthologues-91.

3. Stwórz drzewa filogenetyczne na podstawie macierzy przy pomocy klasy Bio.Phylo.TreeConstruction.DistanceTreeConstructor zarówno metodą UPGMA – hierarchiczną jak i nj (neighbor joining)

4. Wyświetl uzyskane drzewa przy pomocy metody draw_ascii() i draw()

5. Zapisz uzyskane drzewa do formatów newick i phyloxml. obejrzyj wyniki.

W związku z prośbami o plan tematów na kolejnych wykładach i ew. lektury podaję tu plan na dotychczasowe i planowane zajęcia do kolokwium:

Sekwencje DNA i komplementarność, grafy – z książki Pevznera “Computational Molecular Biology”)

Modele ewolucji DNA – część z notatek N. Dojera , część z książki Felsensteina  “inferring phylogenies”)

Porównania par sekwencji – (wg książki Durbina “biological sequence analysis”)

Budowa drzew filogenetycznych – (z książki Felsensteina “inferring phylogenies”)

Uliniowienia wielu sekwencji (z książki Durbin’a – biological sequence analysis)

Ukryte modele Markowa (HMM) – (z książki Durbin’a j.w.)

Zastosowania HMM do modelowania domen białkowych

Wyszukiwanie sekwencji podobnych i algorytm BLAST (część notatek na stronie N. Dojera )

Rodziny i funkcje genów

Uzgadnianie drzew

 

 

Dziś zajmujemy się uliniawianiem par sekwencji – wyk3

Zadania na laboratorium:

1. Zapoznaj się z modułem Bio.pairwise2 i klasą align oraz funkcją format_alignment()
2. Wczytaj sekwencje DNA histonów histones.fa i czynników bZIP bzips.fa do pamięci.
3. Dokonaj porównań pomiędzy sekwencjami DNA białek histonowych i bzip – dla każdej pary policz oceny dla najlepszych globalnych i lokalnych uliniowień z afiniczną funkcją kary. Wylicz średnią ocenę w ramach grupy bzip, w ramach grupy histonów i pomiędzy grupami.
4. Dokonaj tłumaczenia sekwencji DNA na białka, powtórz obliczenia używając macierzy substutucji BLOSUM 50

Praca domowa:

Napisz program, który wylicza optymalne lokalne uliniowienie dla sekwencji DNA dla różnych możliwych tłumaczeń na białka (tzn zakładając standardową tablicę kodonów ) ale możliwe różne ramki odczytu w obu sekwencjach. Dla uproszczenia załóżmy, że insercje i delecje powinny występować tylko “trójkami” nukleotydów. Zastanów się, jakby wyglądał algorytm programowania dynamicznego, gdybysmy rozważali ogólną postać insercji i delecji