APB wykład 4 – debugowanie, profilowanie, logowanie

Dziś mamy wykład zdalny. Spotykamy się o 14:15 (można będzie pewnie dołączać od 14tej) na google meet. Link wyślę do Państwa przez usosmail i rocket chat.

Spróbujemy dziś porozmawiać o debugowaniu, profilowaniu, logowaniu w pythonie na podstawie tych slajdów:wyk4-debugowanie. Jeśli starczy nam czasu, wrócimy też do tematów z poprzedniego wykładu, który się nie odbył, nt. automatycznych testów oprogramowania.

Podczas labów, chciałbym znowu porozmawiać z Państwem nt. postępów w pracy zespołowej. Proponuję, żebyśmy tym razem zrobili to także na platformie google meet, w tym samym miejscu gdzie wykład. Jeśli mają Państwo preferencje co do czasu, proszę zapisywać się (zespołowo) pod tym linkiem: https://doodle.com/poll/pp2nb3vtqnursxzb

Przypominam, że do przyszłej środy spodziewam się od każdego z zespołów opisu Państwa projektu. Zdaję sobie sprawę, że stan zagrożenia epidemiologicznego nie pomaga Państwu w pracy zespołowej, ale jestem przekonany, że mogą Państwo wykorzystać dostępne narzędzia do efektywnego działania. Możemy też poświęcić część wykładu na dyskusję nad problemami technicznymi w komunikacji i pomyśleć nad najlepszymi technologiami, które mogą Państwu pomóc.

WBO – assignment 1 – Coronavirus phylogeny – zadanie 1

As we are spending our time in quarantine, I thought we might try and consider reconstructing some phylogenetic trees that will tell us something about the evolutionary history of the coronavirus 2019-nCov. In particular, we will try and replicate (not completely, but the general idea) the panels b and c from Figure 2 in this paper:Identification_of_a_novel_coronavirus . The idea is to create and compare two phylogenetic trees created from the whole genome sequence of the coronavirus and from its spike protein (you might read more about spike proteins in this article, back from the previous SARS epidemic).

Our goal will be to write a script that creates these two plots using the following steps:

 

  1. Downloading the sequences from Entrez using Bio.Entrez (23 pts)
  2. Creating multiple alignments of both DNA (whole genome) and protein (spike protein) sequences (3 pts)
  3. Creating phylogenetic trees for both protein and genome alignments using the UPGMA and neighbor joining algorithms (3pts)
  4. Visualizing the trees (1 pt)

Now let us describe all these steps in more detail

1. We need to download whole genome sequences and spike protein sequences from the Genbank Entrez databases (“nucleotide” and “protein” databases, respectively), together with some other viruses. We will only look at 6 species:

  • 2019-nCov virus (current coronavirus)
  • SARS – another coronavirus back from 2002
  • bat coronavirus
  • MERS – middle eastern respiratory syndrome coronavirus
  • influenza A –  flu virus, to see how far are coronaviruses from the flu
  • hepatitis A virus – less related virus

You can use the Bio.Entrez interface. An example of how this works can be seen in a jupyter notebook here (please remember to use your e-mail and not mine in your code…). We will need to do this for all 6 species. You should provide a single script that does that for all species, using different query terms for each genome and protein. You will probably need some trial and error work here, but I want to see just the end result – a script that works.

2. Here we will want to run the ClustalW algorithm on the downloaded sequences. You will need to have the clustalw executable installed (it is the package clustalw in ubuntu). Your script should use the Bio.Align.Applications.Clustalw interface. in the result you should get two multiple sequence alignments (for genomes and for spike proteins).

3. Here you should use the two multiple alignments to create distance matrices and from these distance matrices to create phylogenetic trees (both with upgma and nj methods).

4. You should create images of the 4 trees (upgma or nj x protein or genome) using the draw function.

A full solution is a script (in a .py file sent as a real attachment -not some usos attachment, not an ipynb, not a link) together with the images for part 4 (in png format). All of this should be sent to me (bartek[at]mimuw.edu.pl) with the tag [WBO] in the e-mail subject by April 6th. The code should be reasonably easy to read (with comments and the code divided into sensible functions or classes. Not a huge spaghetti of code.)

 

WBO 4 – phylogenetic trees

Today we will hold the lecture remotely via google meet. I will send you the link via the chat server. If you cannot access the chat, please let me know by e-mail.

Polish version below…

Today we will discuss the simple methods for phylogenetic tree reconstruction based on a distance matrix. The slides are in two parts: My part regarding simple progressive heuristic methods: wyk4  and another part on ML methods by P. Górecki : gorecki-ml-tree

During the lab, we will try to build trees from sequences in practice. First, we should get ourselves acquainted with the Bio.Phylo module  tutorial  and wiki.

1. Let us first consider two lists of sequences containing paralogues and orthologues of phenylalanine hydroxylase gene as well as histone H2B paralogues and orthologues. Let us assume, that we have them already aligned (we will cover the multiple alignment process next week). See the files  Human_PAH_orthologues,  Human_PAH_paralogues and  Human_H2BFS_paralogues ,  Human_PAH_orthologues-91. Let  us use the read() function from the Bio.AlignIO module to read these files.

2.Calculate distance matrices for these groups of sequences using the  Bio.Phylo.TreeConstruction.DistanceCalculator  (using  BLOSUM62 substitution matrix, separately for each file).

3. Create the trees using Bio.Phylo.TreeConstruction.DistanceTreeConstructor both for UPGMA and nj (neighbor-joining)

4. look at the resulting trees using  draw_ascii() and draw() methods

5. Write out the resulting trees into newick and phyloxml formats.Compare resulting files.

Wersja polska:

Dziś rozmawiamy o prostych metodach kostrukcji drzew na podstawie macierzy odległości. Slajdy na dziś składają się z mojej prezentacji nt. metod heurystycznych:  wyk4 i slajdów P. Góreckiego nt. metod ML: gorecki-ml-tree

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 listy sekwencji. Rozważmy sekwencje paralogiczne i ortologiczne ludzkiej hydroksylazy fenyloalaniny.  Załóżmy, że mamy już sekwencje aminokwasowe uliniowione globalnie.  Ze względu na czasochłonność procesu uliniowienia, użyjemy plików zuliniowieniami  Human_PAH_orthologues i Human_PAH_paralogues oraz Human_H2BFS_paralogues i (jeśli ktoś liczy lokalnie – uwaga duży) plik Human_PAH_orthologues-91. Wczytaj te pliki przy pomocy metod z modułu Bio.AlignIO

2.Wylicz macierze odległości dla tych grup sekwencji przy pomocy klasy Bio.Phylo.TreeConstruction.DistanceCalculator  (dla macierzy BLOSUM62, osobno dla paralogów i osobno dla ortologów genów PAH – to  zajmie chwilę).

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.

ONA 3 – Wykresy

Dzisiaj zajmujemy się wykresami. Wykład jest w formie google meet, link umieszczę na naszym chacie https://chat.mimuw.edu.pl/channel/obliczenia-naukowe-bibs . Osoby nie posiadające dostępu do chatu, proszę o informację mailem, spróujemy znaleźć inne rozwiązanie. Laboratorium także będzie w formie zdalnej, pytania można zadawać przez chat.

Slajdy są tutaj ONA3-wykresy

Zadania na dziś:

Continue reading “ONA 3 – Wykresy”

ONA 2 – Wektory, macierze i numpy

Slajdy do wykładu – ONA2-Macierze

Jako interpretera pythona mogą Państwo dalej używać naszego serwera Jupyter. Proszę o nie nadużywanie tych interpreterów – szczególnie jeśli chodzi o zajmowanie pamięci. Trochę porad jak używać jupytera znajdą Państwo tu .

Ćwiczenia:

  1. Stwórz macierze A i B o wymiarach  1000×1000 zawierające wartości A[i,j]==i*3-j*5 i B[i,j]==np.sqrt(A[i,j]) (Zauważ, że to wymaga wartości urojonych, a więc macierzy typu complex64 i wartość sqrt(-1) ->.j)
  2. Napisz funkcję poz(indeksy, shape), która dla zadanej krotki indeksów zwraca liniową pozycję elementu o zadanych indeksach w wielowymiarowaj macierzy o kształcie shape. Np. poz([1,1,1],(2,2,2)) daje 7, zaś poz([0,1,0],(2,2,2)) daje 2. Jeśli ktoś przyśle mi w ciągu tygodnia taki program liczący poz() dla dowolnej liczby wymiarów, i nie tworzącej dodatkowych macierzy, e-mailem, z dopiskiem [ONA] w temacie, to otrzyma dodatkowy punkt do zaliczenia.
  3. Napisz program, który  zmienia znak elementów w macierzy B, których suma indeksów jest nieparzysta, na ujemny. Warto wykorzystać  mnożenie macierzy i ew. broadcasting.
  4. korzystając z metody sort() , posortuj macierz 3-wymiarową wg. drugiej współrzędnej i obejrzyj wynik ze zrozumieniem tego co się wydarzyło.
  5. (*)Korzystając z funkcji frompyfunc() napisz funkcje wektorowe, które liczą a) sumę kwadratów dwóch macierzy  (dwie macierze do jednej)
    b) iloraz i resztę z dzielenia całkowitoliczbowego przez 17 (jedna macierz wejściowa i dwie macierze wyjściowe)
  6. Korzystając z modułu time, porównaj prędkość funkcji liczących N potęg dwójki – na liście i w wektorze. Dla jakich wartości N warto używać macierzy?

APB 1,2 – Wprowadzenie i projekty BioX

Na pierwszym wykładzie rozmawialiśmy głównie o planach na wykład. Plan na lab, to zapoznanie się z podstawami systemów kontroli wersji.

Slajdy są tu: wyk1-intro

Materiały do labu:

Na drugim wykładzie rozmawialiśmy o bioperlu, biopythonie i innych projektach open-source. Slajdy do wykładu są tu: wyk2-bioX

Zgrubny harmonogram pracy na laboratoriach (zachęcam do wyprzedzania):

26 II – zebranie zespołów, (ew. wybór kierowników, może być później)
11 III – decyzja dot. tematu projektu
 praca w grupach
25 III – określenie zakresu projektu – krótki opis projektu (dostarczony do mnie dokument 1-2 stronicowy
 praca w grupach
15 IV – szczegółowy opis funkcjonalności do stworzenia przez poszczególne osoby
 praca w grupach
29 IV – prezentacja początkowych implementacji, dyskusja decyzji projektowych
 praca w grupach
4 VI – prezentacje projektów

 

WBO 2 – Sequence evolution models

Today, we will discuss basic mathematical models of sequence evolution. (polska wersja poniżej)

The slides for the lecture can be found here: wyk2

Again, you can use our jupyter server. There are some tutorials on jupyter usage  here. Otherwise, you will need a python interpreter with biopython installed.

Assignments for today:

1. Take a look at the Bio.SubsMat.MatrixInfo module, where PAM and BLOSUM matrices are defined.

2. Implement a simple, discrete Markov model, that allows for simulation of mutations in a sequence over discrete time steps, according to substitution matrices defined using Bio.SubsMat. Try to implement Kimura and Felsenstein models.

3. Take the test_fasta file from last week and see how many steps of random mutations without selection does it take to change one of the sequences to another without selection pressure.

4. Implement a simple Markov process with selection, i.e. after making a random mutation that brings you closer (according to Hamming distance) to the target sequence, then accept it. Reject it otherwise. Compare run times of this implementation with the previous one (without selection).

 

Dziś będziemy mówić o modelach ewolucji sekwencji.

Slajdy do wykładu są tu: wyk2.

Podczas labów mogą Państwo korzystać z serwera Jupyter, który ma zainstalowany pakiet biopython. (kilka porad jak używać jupyter’a tu ). Zadania do wykonania na dziś:

1. Obejrzyj moduł Bio.SubsMat.MatrixInfo, gdzie zdefiniowane są modele PAM i BLOSUM

2. Zaimplementuj prosty łańcuch Markowa, który pozwoli na symulowanie mutacji w czasie dyskretnym, zgodnie z macierzami przejścia zdefiniowanymi wykorzystując moduł SubsMat (substitution table). Rozważ model analogiczny do modelu Kimury i Felsensteina

3. Rozważ sekwencje z pliku fasta z poprzednich zajęć. Zbadaj ile mutacji losowych potrzeba aby Łańcuch Markowa “przeszedł” pomiędzy dwiema sekwencjami bez selekcji.

4. Zaimplementuj prosty proces Markowa z selekcją – jeśli mutacja zbliża nas do “celu” w sensie odległości Hamminga, to ją akceptuj, w przeciwnym razie wracaj do punktu wyjścia.

 

ONA 1 – arytmetyka zmiennoprzecinkowa

Slajdy do wykładu można znaleźć tu: ONA1-Arytmetyka

 

Uwaga nt. liczb całkowitych. Obecnie python3 ukrywa zupełnie przed nami fakt istnienia ograniczeń wielkości rejestrów przy operacjach na liczbach całkowitych. Aby zobaczyć jak to działa, warto uruchomić interpreter python2.7 i wykonać operacje na dużych liczbach.

Pracujemy w narzędziu python3 (nie w octave) i wykonujemy zadania częściowo inspirowane zadaniami z labu 3 do metod numerycznych. Zachęcam do użycia interpretera ipython3. Nieco więcej technicznych szczegółów jest też dostępne w materiałach do wykładu z metod numerycznych.

Będziemy też korzystać z modułu numpy (pisząc “import numpy as np”)  a w szczególności z jego typów danych float64 i float128

Jeśli ktoś ma kłopoty z interpreterem, to może korzystać z serwera jupyter, używając hasła  ‘2020-studenci’. Proszę wtedy pracować w osobnych podpisanych  notebookach.

W szczególności:

1. Policz epsilon maszynowy (zarówno jako wartość dziesiętną ale też liczbe bitów cechy binarnej przy pomocy sprawdzania dla jakich n 1.0+2.0**(-n) >1.0. Co by się stało, gdybyśmy wykonali tę próbę w okolicy 0.0 zamiast 1.0. Czy taki epsilon maszynowy zależy od liczby bitów cechy, czy mantysy?

2. Oblicz wartość funkcji (pochodzącej od Siegfrieda Rump’a): f(x,y)= 333.75*y**6+x**2*(11*x**2*y**2-y**6-121*y**4-2)+5.5*y**8+x/(y*2) w punkcie x=77617, y=33096. Wykonaj obliczenia w pythonie2, pythonie3, np.float64, np.float32, np.float128 i module Fractions. Jakie błędy obliczeń tam zachodzą? Jeśli ktoś napisze program obliczający tę wartość prawidłowo, i prześle mi e-mailem program i prawidłowy wynik, to dostanie dodatkowy punkt do zaliczenia.

3(*). Zastanów się nad zadaniem z szeregami zbieżnymi (lub rozbieżnymi) numerycznie ( w materiałach do labu 3). Zaimplementuj liczenie sum S(n) w pythonie, dla dowolnego n i dla kilku szeregów podanych w rozwiązaniu. Które z nich można dokładnie policzyć pry użyciu modułu Fractions?

4(*). Policz wartości e(x) wg rozwinięcia wzoru Taylora dla zadanego maksymalnego stopnia. Porównaj różnicę wyniku jaki otrzymujesz dla x=-90 i x=2.

5(*). Policz sumy częściowe szeregu harmonicznego przy użyciu liczb zmiennoprzecinkowych i przy użyciu modułu Fractions. Porównaj wyniki.