Month: May 2016

Dzisiejszy wykład porównujący cechy różnych pakietów do obliczeń jest tu:ONA11-Pakiety-oprogramowania

Zadania na dziś są dość proste i dotyczą przede wszystkim obliczeń symbolicznych i pakietu sympy:

1. Przy użyciu pakietu sympy, zbadaj podstawowe możliwości obliczeń symbolicznych: stwórz zmienną sumboliczną x (Symbol(‘x’)), stwórz wyrażenie wielomianowe (np. x**2-2) oraz wyrażenie trygonometryczne (np. sin(x)) i policz dla nich pochodne, całki (oznaczone i nieoznaczone) oraz granice (w 0, z prawej i lewej strony oraz , granice 1/x w nieskonczonosci: sympy.oo,).
2. Zobacz jak wyniki tego typu zadań wyglądają w Wolfram Alpha
3. Zarejestruj się w systemie SageMath i tam spróbuj rozwiązać te same problemy
4. (*) Spróbuj rozwiązać symbolicznie układ równań różniczkowych zwyczajnych z poprzednich zajęć (model Lotki-Volterry albo Zombie INvasion) w jednym z wymienionych wcześniej pakietów

 

Dziś będziemy zajmować się różniczkowaniem i całkowaniem numerycznym. Slajdy są tu: ONA10-calkowanie

Przydatny będzie nam przede wszystkim pakiet scipy.integrate ale będziemy też korzystać z nowych funkcji pakietu matplotlib, takich jak wykresy strzałkowe

0. Na początek proponuje obejrzeć sobie program opisujący model Lotki-Volterry lotka_volterra i spróbować wykonać go po kawałku ze zrozumieniem. Jeśli są pytania, to proszę pytać

1. Naszym właściwym zadaniem jest budowa podobnego modelu dla problemu inwazji Zombie (opisanego w artykule naukowym ). Mamy tu 3 zmienne: ludzi niezarażonych (S), zombie (Z) oraz ofiary konfliktu (R). System jest dynamiczny, a jego zmiany opisują równania różniczkowe:

• dS/dt = P – B*S*Z -D*S,
• dZ/dt = B*S*Z +G*R – A*S*Z,
• dR/dt = D*S + A*S*Z – G*R,

gdzie

• P oznacza przyrost naturalny ludzi,
• B – współczynnik z którym zombie spotykając niezarażonego może go zarazić,
• D – współczynnik śmiertelności ludzi z przyczyn naturalnych,
• G – prawdopodobieństwo, że zombie powstanie z martwej ofiary,
• A – szansa na to, że w pojedynku człowieka z zombie, to zombie zostanie pokonane

Interesują nas następujące pytania:

a) Jak wygląda ewolucja systemu w czasie (w zakresie 0..5000),  jeśli na początku jest tylko 1000 ludzi niezarażonych (nie ma zombie, ani ciał) i mamy następujące prametry:

• P=0.0001
• D=0.0001
• B=0.0095
• G=0.0001
• A=0.0001

Rozwiązanie przedstaw na wykresie 3 zmiennych (S,Z i R) od czasu

b) Czy zagładę populacji można zatrzymać poprzez zwiększanie współczynnika urodzin P? Rozwiązanie uzasadnij. Który ze współczynników musiałby wynosić 0, żeby populacja mogła poradzić sobie z zombie?

c) Zaprezentuj pole wektorowe populacji zombie i ludzi względem parametrów A i B?

d*) Zaprezentuj pole wektorowe 3d dla zmiennych S,Z,R

 

Naszym zadaniem jest stworzenie pakietu umożliwiającego obróbkę macierzy kontaktów chromosomowych. Będzie on obejmował 3 zasadnicze zadania:

1. Wczytywanie listy kontaktów dla wybranego chromosomu do macierzy o zadanej rozdzielczości (5pkt). Mamy dany plik wejściowy, który zawiera 4 kolumny:

chromosom1 pozycja1 chromosom2 pozycja2

Każdy taki wpis (a będą ich miliony) opisuje kontakt pojedyńczej nici DNA pochodzącej z chromosomu1 na pozycji1 z nicią pochodzącą z chromosomu2 na pozycji2. Mając do dyspozycji taki plik oraz podany numer chromosomu C i długość segmentu S powinniśmy skonstruować kwadratową macierz symetryczną M kontaktów (rozmiaru N*N, gdzie N= L/S, gdzie L oznacza długość chromosomu C) pomiędzy segmentami zadanego chromosomu. Np. jeśli chromosom ma długość L=1000000 par zasad, a S=1000, otrzymamy macierz M rozmiaru 1000*1000, a wartość M[i,j] będzie oznaczać ile wierszy opisywało kontakt pomiędzy fragmentem [1000*i:1000*(i+1)] a fragmentem [1000*j:1000*(j+1)].

2. Normalizacja macierzy (10pkt). Mając daną macierz kwadratową M o zadanym rozmiarze N*N, często potrzebujemy dokonać jej normalizacji. W tym przypadku interesuje nas rozwiązanie równania M[i,j]= B[i]*[B[j]*T[i,j]+epsilon, gdzie wektor B (rozmiaru N) opisuje multiplikatywne obciążenia poszczególnych segmentów, zaś od Macierzy T (rozmiaru N*N) oczekujemy, że będzie symetryczna i znormalizowana (sum(T[:,i])=1.0 dla każdego i). Daje nam to układ równać, który wymaga rozwiązania iteracyjnego, minimalizującego błąd (epsilon). Wynikiem tej procedury powinna być macierz T, wektor B i wartość epsilon, a dozwolonym parametrem powinna być maksymalna liczba iteracji i maksymalny Epsilon.

3. Wyświetlanie  macierzy kontaktów (5 pkt). Kiedy już mamy macierz kontaktów (M lub T) to chcielibyśmy móc obejrzeć ją na ekranie lub zapisać do pliku. Potrzebny nam będzie program, który na podstawie macierzy i dodatkowych parametrów wyświetla zadaną macierz metodą imshow(). Program powinien pozwalać na następujące opcje:

• określanie podmacierzy do wyświetlenia (tylko prostokąt o współrzędnych x_1:x_2,y_1:y_2)
• Wygładzanie mapy przy pomocy filtra Gaussowskiego o określonym promieniu
• logarytmiczna lub liniowa mapa kolorów
• wybór schematu kolorów
• Zmniejszanie rozdzielczości dużych map do wizualizacji o stały czynnik (jeśli mapa ma np. 10000×10000 można zadać zmniejszenie jej przed wizualizacją – np. zmniejszenie o czynnik 10 powinno uśrednić dane z kwadratów rozmiaru 10×10 i wyświetlić mapę rozmiaru 1000×1000)

Wynik wysyłamy w postaci pojedyńczego skryptu (jako niespakowany załącznik) na adres prowadzącego zajęcia z dopiskiem ONA-3-2016 w tytule do dnia 20. VI 2016.

UWAGA: Przykłądowy plik z danymi (dla 1 chromosomu o długości 1mln par zasad) można znaleźć tutaj – uwaga, plik jest duży…)

UWAGA 2: Kontakty chromosomowe są dość często wzbogacone w kontakty blisko “przekątnej”, czyli kontakty niedużego zasięgu. W związku z tym warto rozważyć wersję programu, która ignoruje wartości w macierzy na samej przekątnej.

Wskazówki do rozwiązania po pierwszych pytaniach:

1) Czy dla dużych macierzy M (np. 1000×1000) możemy myśleć o rozwiązaniu dokłądnym tego układu równań?

2) Czy normalizacja kolumnowo/wierszowa macierzy T nasuwa nam jakieś oszacowanie na przybliżone wartości B[i]?

3) Czy jeśli ktoś dałby nam wartość przybliżoną B, to umielibyśmy znaleźć krok iteracji poprawiający to oszacowanie przy pomocy macierzy M?

4) Epsilon, o którym wspominam w zadaniu jest pewnego rodzaju skrótem myślowym. Kiedy występuje w równaniu, mamy na myśli różnicę między prawą, a lewą stroną tego konkretnego równania (powiedzmy, mały epsilon). Kiedy pojawia się jako kryterium zakończenia iteracji (duży Epsilon), jest to już różnica między wszystkimi lewymi i prawymi stronami równań. W zagadnieniu najmniejszych kwadratów, typowo ten duży Epsilon to suma kwadratów małych epsilonów. W naszym przypadku, gdzie liczba równań, a więc i składników sumy zależy istotnie od parametrów zadanych przez użytkownika, warto rozważyć duży Epsilon będący maksimum z wartości bezwzględnych małych epsilonów.

Zadanie daje Państwu maksymalnie 10 punktów. Realizowane jest w zespołach 2,3, lub 4 osobowych.

Punkty będą przydzielone po prezentacji projektów na ostatnich zajęciach.

Zadanie polega na zaprogramowaniu prostej gry, która wykorzystywać będzie czujnik kinect

Przykładowa gra wygląda tak:

Dwóch użytkowników staje w zasięgu widzenia kamery kinect.

Komputer losuje dwie liczby x_1, x_2z zakresu [0,2000] (tak naprawdę sensownego mniejszego zakresu), i pozwala graczom w turach zgadywać je poprzez ustawianie się w określonej odległości od kamery.

Po każdej turze program wykrywa obecność graczy, mierzy ich średnie odległości od kamery d_1, d_2 i jeśli |x_i -d_i| < epsilon, to gracz i wygrał, a jeśli nie, to komputer wyświetla informację: za blisko, albo za daleko. Jeśli któryś z graczy osiąga “dobrą” odległość wcześniej, to wygrywa, jeśli oboje “trafiają” w tej samej turze, to jest remis.

Oczywiście zachęcam Państwa do wyboru innych ciekawszych dla Państwa scenariuszy gry, ale proszę o zgłaszanie mi ich do akceptacji.

Proszę o podzielenie się na zespoły i wymyślenie projektów do następnych zajęć: 23 maja. Prezentacje będą na ostatnich zajęciach: 13. Czerwca (było omyłkowo 6.). Prezentacje przesunęliśmy na 21. czerwca o 12:00.  Odbywają się w sali 3041 (tam gdzie odbywały się laboratoria).

Wolałbym, aby Państwa programy koncentrowały się na kreatywnej obróbce danych (np. wykrywanie dłoni, czy mierzenie “plam” sylwetki) niż na bardzo szybkiej “zręcznościowej” rozgrywce, ale nie ma formalnych ograniczeń.

 

 

 

Dzisiejsze zajęcia poświęcimy na analizę obrazów 2d. Slajdy są tu ONA9-Obrazy

Przydać mogą się pakiety: scipy.ndimage i matplotlib

Zadania na dziś:

0. Pozyskaj przy pomocy kinect swoje zdjęcie w pracowni komputerowej zarówno w postaci macierzy głebokości, jak i obrazów rgb

1. Przekształć obraz RGB do skali szarości poprzez uśrednienie składowych.

2. Wyświetl te obrazy na swoim komputerze przy pomocy biblioteki matplotlib

3. Wyświetl histogram obrazu i dokonaj jego “wyrównania” pisząc program w  języku python

f(x) = (x-min)/(max-min)*255 dla parametrów min i max

4. Zastosuj filtr Gaussowski rozmiaru k, aby “wygładzić” obraz i obejrzyj wyniki dla różnych k (3,5,7,9,11,…)

5. Napisz program wykrywający krawędzie przy pomocy filtru Sobel’a

6(*). Napisz program wykrywający cienie na obrazie z miernika odległości.

Zaimplementuj procedurę MCL (Markov Chain Clustering) dla zadanej macierzy kwadratowej.

Program powinien być w postaci skryptu, który będzie pozwalał podać nazwę pliku i sposób reprezentacji macierzy M (tekstowy, binarny, macierz rzadka w pliku .mat) oraz parametry r, N i epsilon). Następnie powinien wykonywać iterację wg następującego planu:

• dopóki i<N oraz |M(i)-M(i+2)|>epsilon:
•        M(i+1) = M(i)*M(i) #mnożenie macierzy
•        M(i+2)= Inflate(M(i+1),r)

gdzie funkcja inflate wykonuje operację M[i,j]=M[i,j]**r (po elementach), a następnie normuje kolumny do 1.

Skrypt powinien wypisywać macierz w wybranym formacie (rzadkim, binarnym lub tekstowym) oraz na życzenie wykreślać wykres |M(i+2)-M(i)| względem i.

Warto dodać do programu funkcję, która zwraca spójne składowe grafu opisanego taką macierzą, czyli listy wierzchołków, pomiędzy którymi pozostają połączenia o niezerowej wadze.

Termin oddania prac – 30. maja, prace wysyłamy e-mailem do prowadzącego z dopiskiem ONA-2-2016 w tytule, w postaci pojedyńczego pliku pythona jako załącznik.

Przykładowe dane – macierz podobieństwa szczepów konopii (na podst. pracy  Hakki et al, Electronic Journal of Biotechnology, Vol 10, No 4 (2007)) : cannabis4.npy, Dla tych danych sensownymi wartościami parametrów są N=20, epsilon = 0.5, r=1.1

 

Slajdy do wykładu są tutaj: ONA8-DSP

Interesują nas przede wszystkim biblioteki scipy.signal i scipy.fftpack

1. Wygeneruj sygnał funkcji wielomianowej o zaburzeniu Gaussowskim, dla 10000 punktów. Wykonaj uśrednienia tego sygnału przy pomocy średniej kroczącej, albo splotu przy pomocy filtru kwadratowego, trójkątnego lub Gaussowskiego. WYniki zaprezentuj na wykresie

2. Pobierz sygnał o zajętości nukleosomowej z plików pochodzących z eksperymentu MNase-Seq (tutaj). Wczytaj go do numpy jako wektor. Przedstaw go na wykresie. Używając szybkiej trasformaty fouriera oblicz widmo Fourierowskie tego sygnału. Przedstaw je na innym wykresie. Wyzeruj część widma wysokich częstotliwości i dokonaj odwrotnej transformaty Fouriera, aby uzyskać wygładzony sygnał zajętości nukleosomami. Spróbuj tak dobrać parametry filtra, aby wykres po odwrotnej transformacie Fouriera miał okres zbliżony do 160-200 par zasad.