bioinformatics-basics/laboratorium_3.qmd at main · gladkia/bioinformatics-basics · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
---
title: "Laboratorium 3: Profile sekwencyjne, Ewolucja i ukryte wiadomości"
format: html
engine: knitr
filters:
  - webr
---

```{r, echo = FALSE, results='asis'}
library(webexercises)
```

## Część 1: BLAST, E-value i... Dinozaury!

W 1990 roku Michael Crichton wydał słynną powieść "Park Jurajski". Na stronie 103 przedstawił sekwencję nukleotydową, która rzekomo miała pochodzić z wyekstrahowanego DNA dinozaura.

W prawdziwym świecie naukowcy, tacy jak Mark Boguski, szybko zauważyli pewne nieścisłości. Boguski zaproponował własną, poprawioną sekwencję, w której... ukrył wiadomość! Aby ją odczytać, należy przetłumaczyć sekwencję na białko i przyrównać ją do bazy danych *nr* używając odpowiedniego narzędzia BLAST.

### Jak ocenić wiarygodność znaleziska? (E-value)

Zanim pobawimy się w detektywów, musimy zrozumieć, jak algorytm BLAST ocenia, czy podobieństwo dwóch sekwencji nie jest dziełem przypadku. Służy do tego parametr **E-value** (Expect value):

* Im E-value jest **mniejsze** (bliższe zeru), tym wynik jest **bardziej wiarygodny**.
* Przykładowo, E-value rzędu $6e-145$ uznaje się za ekstremalnie znaczące, podczas gdy $0.39$ sugeruje, że dopasowanie może być przypadkowe.

**Sprawdź swoją wiedzę:** Który z poniższych wyników BLAST jest najbardziej znaczący statystycznie?
`r mcq(c("E-VALUE: 0.39", "E-VALUE: 2e-29", answer = "E-VALUE: 6e-145"))`

---

## Część 2: Multiple Sequence Alignment (MSA) i WebLogo w R

Kiedy mamy do czynienia z wieloma sekwencjami, wykonujemy przyrównanie wielu sekwencji (MSA). Dzięki temu możemy znaleźć regiony konserwatywne – części sekwencji, które nie zmieniają się w toku ewolucji, co sugeruje ich kluczową funkcję.

Graficzną reprezentacją tej konserwatywności jest **Sequence Logo**. Wysokość każdej litery odzwierciedla, jak bardzo dany aminokwas jest konserwatywny na określonej pozycji.

### Interaktywne zadanie: Stwórz własne Logo w R

Użyjemy pakietu `ggseqlogo`, aby zwizualizować zestaw krótkich sekwencji peptydowych.

**Zadanie:** Uruchom poniższy kod. Zwróć uwagę, na których pozycjach aminokwasy są wysoce konserwatywne, a na których ulegają mutacjom.

```{webr-r}
# Instalacja i ładowanie bibliotek do wizualizacji
webr::install("ggplot2")
webr::install("ggseqlogo")
library(ggplot2)
library(ggseqlogo)

# Nasze przyrównanie wielu sekwencji (MSA) w formie wektora
alignment <- c(
  "VDFWAE",
  "VDFWAP",
  "VDFWAE",
  "VDFWAE",
  "VDFWAP",
  "VDFSAT",
  "VDFSAT",
  "VDFYAT"
)

# Rysujemy piękne logo!
ggseqlogo(alignment) +
  ggtitle("Sequence Logo naszych peptydów") +
  theme_minimal()
```

---

## Część 3: PSI-BLAST i Profile PSSM

Zwykły BLAST może przeoczyć dalekich homologów. Wtedy wkracza **PSI-BLAST** (Position-Specific Iterative BLAST). Jak to działa?

1. Wykonuje wstępne przeszukiwanie.
2. Z najlepszych wyników buduje profil **PSSM** (Position-Specific Scoring Matrix).
3. Używa tego profilu, by szukać głębiej, powtarzając proces (iteracje) aż do osiągnięcia **konwergencji** (gdy przestają pojawiać się nowe sekwencje).

**Zadanie końcowe (Detektywistyczne):**
Wykonaj kilka iteracji narzędzia PSI-BLAST na stronie NCBI używając sekwencji Marka Boguskiego (plik `dino2.txt`).

* Sprawdź wyniki przyrównania, aby znaleźć ukrytą wiadomość.
* Do jakiej rodziny białek tak naprawdę należy ta sekwencja?