Motywy i pętle
Strona 1
Motywy i pętle
Bardzo często zdarza się, że chcemy odnaleźć w sekwencji DNA lub białka jakiś interesujący nas motyw, czyli inaczej
mówiąc jakiś określony fragment sekwencji. Perl oczywiście daje nam narzędzia do zrealizowania tych celów.
Zaczniemy jednak od…
Wyrażenia warunkowe
Jeżeli ktoś miał styczność z językiem BASIC z pewnością nie jest mu obca struktura IF-THEN. Wyrażenia warunkowe
pozwalają na dokonywanie wyborów zależnie od wyników przeprowadzonego testu. W Perlu występuje kilka rodzajów
wyrażeń warunkowych. Są to if, if-else i unless. W tym momencie warto zwrócić uwagę na jeden szczegół. W
konstrukcjach warunkowych Perla do przeprowadzenia testów prawda-fałsz, czyli mówiąc krótko do sprawdzenia, czy
spełniony jest interesujący nas warunek, używa się operatora „==”. Użycie pojedynczego znaku „=” byłoby błędem,
ponieważ „=” jest zarezerwowane do operacji przypisywania wartości zmiennym.
Wyrażenie if:
if(1 == 1) {
print „1 równe jest 1\n\n”;
}
daje w wyniku wydruk na ekranie „1 równe jest 1”, ponieważ w rzeczy samej jedynka równa jest jedynce.
Wyrażenie:
if(0 == 1) {
print „0 równe jest 1\n\n”;
}
nie da nam żadnego wyniku, ponieważ warunek nie jest spełniony. Zero nie równa się jedynce.
Wyrażenie tego samego typu można zapisać w inny sposób: print „1 równe jest 1” if (1 == 1); Wynik
będzie taki sam jak w pierwszym przypadku.
Wyrażenie if-else pozwala na wykonanie dwóch różnych działań w zależności od wyniku testu.
if(1 == 0) {
print „1 równe jest 0\n\n”;
} else {
print „1 nie jest równe 0\n\n”;
}
Po uruchomieniu skryptu otrzymamy wydruk „1 nie równa się 0”, zgodnie z prawdą.
Wyrażenie unless oznacza dokładnie to samo, co w języku angielskim, a mianowicie „if not”.
unless(1 == 0) {
print „1 nie równa się 0\n\n”;
}
W wyniku otrzymujemy wydruk „1 nie równa się 0”, ponieważ w rzeczy samej test warunkowy zwrócił nam prawdę, 1
nie jest równe 0.
Warto zwrócić uwagę na symbole { }. Pomiędzy nimi znajduje się część kodu, która jest wykonywana w przypadku, gdy
test warunkowy zwrócił nam w wyniku prawdę.
Motywy i pętle
Strona 2
Sprawdzenie równości dwóch liczb nie jest jedynym warunkiem, który możemy podać w teście warunkowym. Inne
opcje to sprawdzenie, czy dane elementu nie są sobie równe z użyciem symbolu „!=”, czy jeden jest mniejszy od
drugiego „<”, albo większy „>”. W celu porównania łańcuchów tekstowych lub zmiennych tekstowych stosuje się
operatory eq (odpowiednik ==) oraz ne (odpowiednik !=).
Przykład 1. if-elsif-else
#!/usr/bin/perl –w
# if-elsif-else
$word = ‘MNIDDKL’;
if ($word eq ‘QSTVSGE’) {
print „QSTVSGE\n”;
}
elsif ($word eq ‘MRQQDMISHDEL’) {
print „MRQQDMISHDEL\n”;
}
else {
print „Czy \”$word\” to w ogóle peptyd?\n”;
}
exit;
Zwróćmy uwagę na symbol \” w bloku else. Postawienie \ przed znakiem cudzysłowia informuje komputer, że należy
wydrukować ten znak na ekranie, a nie jest to zakończenie tekstu drukowanego przez polecenie print.
Pętle
Pętle pozwalają na powtarzające się wykonywanie bloku poleceń ograniczonego symbolami { }. W Perlu występuje
kilka rodzajów pętli: while, for, foreach i inne. Przedstawiony poniżej przykład ilustruje zastosowanie pętli
while.
Przykład 4. Odczytywanie sekwencji z pliku, podejście 4.
#!/urs/bin/perl –w
# Odczytywanie sekwencji z pliku, podejście 4
$proteinfilename = ‘bialko.pep’;
# Najpierw otwieramy plik do odczytu, a jeżeli to się nie powiedzie komputer
# drukuje komunikat o błędzie i kończy program
unless(open(PROTEINFILE, $proteinfilename)) {
print „Nie mogę otworzyć pliku $proteinfilename!\n”;
exit;
}
# Odczytujemy sekwencję z pliku korzystając z pętli while
# drukujemy każdą odczytaną linie
while($protein = <PROTEINFILE>) {
print „### Oto następna linia z pliku: \n”;
print $protein;
}
Motywy i pętle
Strona 3
close PROTEINFILE;
exit;
W pętli while odczytywana jest kolejna linia z pliku i przypisywana jest ona zmiennej $protein. Na początku pętli
sprawdzane jest, czy można odczytać linie i jeżeli jest to możliwe, to znaczy test warunkowy zwróci w wyniku prawdę,
to wykonywane są instrukcje znajdujące się pomiędzy symbolami { } (ciało pętli). Warto również zwrócić uwagę na
wykorzystanie instrukcji warunkowej unless. W przypadku, gdy nie można otworzyć pliku drukowany jest komunikat
o błędzie, a program kończy się.
Układ kodu
Wraz z rozpoczęciem używania pętli i instrukcji warunkowych warto zastanowić się nad formatowaniem kodu pisanego
przez nas programu. W gruncie rzeczy komputerowi nie zależy na tym, żeby ów kod był przejrzysty i zrozumiały. Ważne,
żeby był prawidłowy. Tymczasem dla nas przejrzystość jest ważna. Dopóki nasze programy zawierają 10 – 15 linii
problemu większego nie ma. Gorzej zaczyna się robić, kiedy ilość linii kodu przyrasta, a pętle i instrukcje warunkowe
zaczynają się zagnieżdżać…
Oto przykłady formatowania kodu z instrukcją if w pętli while. Wszystkie są poprawne i działają.
Format A
while ($zywy) {
if ($głodny) {
print „Wymagany jest pokarm\n”;
}
}
Format B
while ($zywy)
{
if ($glodny)
{
print „Wymagany jest pokarm\n”;
}
}
Format C
while ($zywy)
{
if ($glodny)
{
print „Wymagany jest pokarm\n”;
}
}
Format D
while($zywy) {if($glodny) {print „Wymagany jest pokarm\n”;}}
Formatowanie kodu jest oczywiście kwestią osobistych preferencji, jednak sposoby A i B są polecane bardziej niż C, a
tym bardziej niż D.
Motywy i pętle
Strona 4
Znajdowanie motywów
Jedną z częściej wykonywanych operacji w bioinformactyce jest znajdywanie motywów, krótkich segmentów DNA lub
białek, które nas interesują. Perl posiada zestaw poręcznych narzędzi do odnajdywania motywów w łańcuchach
tekstów, co zaraz zilustrujemy przykładami.
Przykład 3. Poszukiwanie motywów
#!/usr/bin/perl –w
# Poszukiwanie motywów
# Program zapyta użytkownika o nazwę pliku zawierającego sekwencję białka
print „Proszę podać nazwę pliku zawierającego sekwencję białka: „;
$proteinfilename = <STDIN>;
unless(open(PROTEINFILE, $proteinfilename)) {
print „Nie mogę otworzyć pliku \”$proteinfilename\”\n\n”;
exit;
}
# Wczytujemy sekwencję do tablicy
@protein = <PROTEINFILE>;
close PROTEINFILE;
# Zamieniamy tablicę na jedną zmienną tekstową
# Łatwiej jest szukać motywów w jednym łańcuchu, a nie w tablicy łańcuchów
$protein = join(‘’, @protein);
# Usuwamy spacje
$protein =~ s/\s//g;
# Korzystamy z pętli i poszukujemy motywów, które będą wprowadzone z klawiatury.
# Gdy użytkownik nie wprowadzi żadnego motywu program się kończy.
do {
print „Wprowadź motyw, którego chcesz szukać: „;
$motif = <STDIN>;
# Usuwamy znak nowej linii z motywu
chomp $motif;
# Szukamy motywu
if ($protein =~ /$motif/) {
print „Znalazłem!\n\n”;
} else {
print „Nie znalazłem\n\n”;
}
# Program zakończy się, gdy użytkownik nie poda motywu
} until ($motif =~ /^\s*$/);
exit;
W programie przedstawionym powyżej znalazło się kilka nowych elementów.
Motywy i pętle
Strona 5
Pobieranie danych od użytkownika: $proteinfilename = <STDIN>. STDIN jest skrótem od standard input i jest
to specjalny rodzaj uchwytu pliku, który jest związany z klawiaturą i umożliwia odczytywanie danych, które zostały z
niej wprowadzone.
Zamiana tablicy na skalar używając join. W linii $protein = join(‘’, $@protein) połączyliśmy ze sobą
wszystkie elementy tablicy @protein i zapisaliśmy je w zmiennej $protein. Poszczególne elementy składowe
nowego łańcucha rozdzielone są znakiem, który umieściliśmy pomiędzy dwa apostrofami (dwa cudzysłowie też by
działały). W naszym przypadku nie było tam nic, więc w łańcuchy zostały połączone bez żadnego specjalnego znaku
lub tekstu rozdzielającego.
Pętla do-until. Pętla ta wykonuje kod zawarty pomiędzy znakami { } dopóki spełniony jest warunek postawiony w
linii until.
Wyrażenia regularne. Pozwalają one na łatwe manipulowanie łańcuchami wszelkich rodzajów. Wyrażenie regularne
jest tak jakby uniwersalną postacią poszukiwanego przez nas motywu, umożliwiającą dopasowanie jednego lub
większej ilości łańcuchów przy wykorzystaniu do ich zapisu specjalnych symboli. Mogą być one bardzo proste, jak całe
słowo, np. /bioinformatyka/ lub bardziej skomplikowane.
W linii $protein =~ s/\s//g; wykorzystaliśmy właśnie wyrażenie regularne. Znany już operator wiązania „=~”
pozwala na przeprowadzenie operacji podstawienia w łańcuchu przechowywanym w zmiennej $protein. Po prawej
stronie operatora pojawiło się wyrażenie regularne „\s”. Symbol ten jest jednym z wielu tak zwanych meta znaków i
oznacza każdy pusty znak, to znaczy spację, znak tabulacji, znak nowej linii itd. W naszym przypadku puste znaki
zastępujemy brakiem znaku, to znaczy po prostu je usuwamy.
Linia if ($motif =~ /^\s*$/) { zawiera kolejne wyrażenie regularne. Tym razem dzięki operatorowi wiązania
możemy sprawdzić, czy wyrażenie to jest obecne w łańcuchu zapisanym w zmiennej $motif. Samo wyrażenie
regularne to: /^\s*$/. Przetłumaczyć je można jako: „dopasuj łańcuch, który od początku (oznaczonego przez ^)
zawiera zero lub więcej (oznaczone *) pustych znaków (oznaczone przez \s) aż do końca (oznaczone przez $)”.
Ostatnie wyrażenie regularne znajduje się w linii if ($protein =~ /$motif/) {. Jak widać komputer znowu
poszukuje motywów w łańcuchu przechowywanym w zmiennej $protein. Tym razem wyrażeniem regularnym jest
tekst, który wprowadził użytkownik i zachowany on został w zmiennej $motif.
Liczenie nukleotydów
Teraz spróbujemy napisać program, przy pomocy którego będziemy mogli policzyć ilość poszczególnych nukleotydów
w sekwencji DNA. Taka informacja może być przydatna na przykład do obliczania procentu par GC w danym genie. Oto
pseudokod naszego programu.
wczytaj sekwencję DNA z pliku
połącz linie w jeden ciąg $DNA
# utwórz tablicę ze zmiennej $DNA
@DNA = explode $DNA
# zainicjuj liczniki
count_of_A = 0
count_of_C = 0
count_of_G = 0
count_of_T = 0
dla każdej zasady w @DNA
Motywy i pętle
Strona 6
jeżeli zasada to A
count_of_A = count_of_A +1
jeżeli zasada to A
count_of_C = count_of_C +1
jeżeli zasada to A
count_of_G = count_of_G +1
jeżeli zasada to A
count_of_T = count_of_T +1
zrobione
wydrukuj count_of_A, count_of_C, count_of_G, count_of_T
Pseudokod pozwala na ułożenie i przeanalizowanie algorytmu, który planujemy zastosować do rozwiązywania naszego
problemu, bez martwienia się o składnię wpisywanych poleceń. Krótko mówiąc jest to plan naszego programu.
Przykład 4. Określanie częstości występowania nukleotydów
#!/usr/bin/perl –w
# Określanie częstości występowania nukleotydów
print „Podaj nazwę pliku zawierającego sekwencję DNA: „;
$dna_filename = <STDIN>;
chomp $dna_filename;
unless(open(DNAFILE, $dna_filename)) {
print „Nie mogę otworzyć plku \”$dna_filename\”\n\n”;
exit;
}
@DNA = <DNAFILE>;
close DNAFILE;
$DNA = join(‘’,@DNA);
$DNA =~ s/\s//g;
# Teraz zamieniamy sekwencję zachowaną w zmiennej $DNA na tablicę
# tak, że każdy element tablicy zawiera jedną zasadę
@DNA = split(‘’,$DNA);
# inicjujemy liczniki
$count_of_A = 0;
$count_of_C = 0;
$count_of_G = 0;
$count_of_T = 0;
$errors = 0;
# Używając pętli liczymy ilość poszczególnych zasad
foreach $base (@DNA) {
if ($base eq ‘A’) {
++$count_of_A;
}
elsif ($base eq ‘C’) {
Motywy i pętle
Strona 7
++$count_of_C;
}
elsif ($base eq ‘G’) {
++$count_of_G;
}
elsif ($base eq ‘T’) {
++$count_of_T;
} else {
print „Błąd, nie znam zasady: $base\n”;
++$errors;
}
}
# Drukujemy wyniki
print „A = $count_of_A\n”;
print „C = $count_of_C\n”;
print „G = $count_of_G\n”;
print „T = $count_of_T\n”;
print „Błędy = $errors\n”;
exit;
Zwróćmy uwagę na nowe szczegóły. Linia @DNA = split(‘’,$DNA); zamienia zmienną $DNA na tablicę @DNA.
Funkcja split jest przeciwieństwem join i również w niej występują apostrofy. Pomiędzy apostrofami możemy
umieścić znak rozdzielający, po napotkaniu którego funkcja będzie tworzyć kolejne elementy. Wstawienie ‘’, czyli
brak znaku rozdzielającego oznacza, że elementy tablicy będą zawierać pojedyncze znaki. W powyższym przykładzie
pojawiła się również inicjalizacja zmiennych, to znaczy tworzenie zmiennej i przypisywanie jej wartości, np.
$count_of_A = 0.
W przykładzie pojawiła się również nowa pętla , foreach. Jest to bardzo wygodna konstrukcja, pozwalająca na
powtarzanie bloku kodu dla wszystkich elementów występujących w tablicy, bez potrzeby znania ich ilości.
W przykładzie 4. znalazło się również wyrażenie typu ++$count. Jest to wygodny sposób zwiększania wartości
zmiennej $count i odpowiada ono wyrażeniu $count = $count + 1, lub w wersji skróconej $count += 1.
Słów kilka o liczbach. Perl w dosyć sprytny sposób podchodzi do zmiennych skalarnych. Przykładowo zmienną można
zapisać jako liczbę l1234 lub jako łańcuch ‘1234’. Perl traktuje obie te zmienne jako łańcuchy podczas drukowania i
jako liczby podczas wykonywania działań arytmetycznych. Zilustruje to następny przykład.
Przykład 5. Liczby i ciągi w Perlu
#!/usr/bin/perl –w
# Liczby i ciągi w Perlu
$liczba = 1234;
$ciag = ‘1234’;
# drukujemy zmienne
print $liczba, „ „, $ciag, „\n”;
# dodajemy zmienne jako liczby
$liczba_lub_ciag = $liczba + $ciag;
Motywy i pętle
Strona 8
print $liczba_lub_ciag, „\n”;
# połącz zmienne jako ciągi
$liczba_lub_ciag = $liczba . $ciag;
print $liczba_lub_ciag, „\n”;
exit;
Operacje na łańcuchach
Aby dostać się do poszczególnych znaków łańcucha nie trzeba go od razu zamieniać na tablicę w sposób pokazany w
przykładzie 4. Co więcej, takiego rozwiązania należy unikać, ponieważ zabiera ono dużo pamięci. Jeżeli tekst jest krótki
to nie ma problemu, jednak jeżeli chcielibyśmy zrobić tablicę zawierającą na przykład genom komara to mogłoby być
trochę gorzej.
Przykład 6. Określanie częstości występowania nukleotydów, podejście 2
#!/usr/bin/perl –w
# Określanie częstości występowania nukleotydów
print „Podaj nazwę pliku zawierającego sekwencję DNA: „;
$dna_file = <STDIN>;
chomp $dna_file;
# Czy plik istnieje?
unless (-e $dna_file) {
print „Plik \”$dna_file\” nie istnieje.\n”;
exit;
}
# Czy można otworzyć plik?
unless (open(DNAFILE, $dna_file)) {
print „Nie mogę otworzyć pliku \”$dna_file\”.\n”;
exit;
}
@DNA = <DNAFILE>;
close DNAFILE;
$DNA = join(‘’, @DNA);
$DNA =~ s/\s//g;
$count_of_A = 0;
$count_of_C = 0;
$count_of_G = 0;
$count_of_T= 0;
$errors = 0;
# Tworzymy pętlę, która sprawdza kolejne nukleotydy w ciągu DNA
for ($position = 0; $position < length $DNA; ++$position) {
$base = substr($DNA, $position, 1);
if ($base eq ‘A’) {
Motywy i pętle
Strona 9
++$count_of_A;
} elsif ($base eq ‘C’) {
++$count_of_C;
} elsif ($base eq ‘G’) {
++$count_of_G;
} elsif ($base eq ‘T’) {
++$count_of_T;
} else {
print „Błąd, nie znam zasady: $base\n”;
}
}
# Drukujemy wyniki
print „A = $count_of_A\n”;
print „C = $count_of_C\n”;
print „G = $count_of_G\n”;
print „T = $count_of_T\n”;
print „Błędy = $errors\n”;
exit;
W kodzie pojawiła się linia unless (-e $dna_filename) {, która sprawdza, czy plik o podanej nazwie istnieje
(operator „-e”). Zauważyć także można kolejną pętlę, for. Blok kodu
for ($position = 0; $position < length $DNA, ++$position) {
# tutaj zrób coś
}
odpowiada pętli while:
while ($position < length $DNA) {
# tutaj zrób coś
++$position;
}
Funkcja length zwraca długość łańcucha, który podaliśmy jej jako argument.
Funkcja substr w linii $base = substr ($DNA, $position, 1); wybiera z ciągu przechowywanego w
zmiennej $DNA fragment o długości 1 znaku rozpoczynając od pozycji $position.
Zapisywanie do plików
Jeżeli potrafimy odczytać zawartość pliku warto byłoby również móc zapisać dane do pliku. Perl oczywiście daje taką
możliwość co za chwilę zobaczymy w kolejnym przykładzie.
Przykład 7. Określenie częstości występowania nukleotydów, podejście 3
#!/usr/bin/perl –w
# Określanie częstości występowania nukleotydów
print „Podaj nazwę pliku zawierającego sekwencję DNA: „;
$dna_filename = <STDIN>;
chomp $dna_filename;
# Czy plik istnieje?
unless (-e $dna_filename) {
Motywy i pętle
Strona 10
print „Plik \”$dna_filename\” nie istnieje.\n”;
exit;
}
# Czy można otworzyć plik?
unless (open(DNAFILE, $dna_filename)) {
print „Nie mogę otworzyć pliku \”$dna_filename\”.\n”;
exit;
}
@DNA = <DNAFILE>;
close DNAFILE;
$DNA = join(‘’, @DNA);
$DNA =~ s/\s//g;
# Inicjujemy liczniki
$a = 0; $c = 0; $g = 0; $t = 0; $e = 0;
# Używamy wyrażenia regularnego i pięciu pętli while
# do obliczania ilości odpowiednich zasad i błędów
while($DNA =~ /a/ig) {$a++};
while($DNA =~ /c/ig) {$c++};
while($DNA =~ /g/ig) {$g++};
while($DNA =~ /t/ig) {$t++};
while($DNA =~ /[^acgt]/ig) {$e++};
print „A=$a, C=$c, G=$g, T=$t, błędy=$e\n”;
# Zapisujemy wynik w pliku „ilość_zasad”
$outputfile = „ilość_zasad.txt”;
unless(open(COUNTBASE, „>$outputfile”)) {
print „Nie mogę otworzyć pliku \”$outputfile\” do zapisu\n\n”;
exit;
}
print COUNTBASE „A=$a, C=$c, G=$g, T=$t, błędy=$e\n”;
close COUNTBASE;
exit;
Pętla while($dna =~ /a/ig) {$a++} zawiera test warunkowy w którym występuje wyrażenie regularne /a/ z
modyfikatorami ig. Modyfikator „i” powoduje, że dopasowywane są zarówno litery małe jak i wielkie, modyfikator
„g” powoduje przeszukiwanie całego łańcucha. Blok {$a++} jest wykonywany w przypadku gdy test warunkowy
zwrócił prawdę i po prostu podnosi o jeden wartość zmiennej $a.
Polecenie open(COUNTBASE, „>$outputfile”) otwiera do zapisu plik, którego nazwa przechowywana jest
w podanej zmiennej. Sam zapis odbywa się za pomocą polecenia print, z tym, że w tym przypadku precyzujemy, że
drukowanie ma się odbywać do uchwytu pliku, czyli krótko mówiąc dane mają być zapisane w pliku.
Motywy i pętle
Strona 11
Obliczenia ilości zasad poszczególnych rodzajów można dokonać również za pomocą operatora transliteracji tr.
Dokonuje on nie tylko podstawienia odpowiednich liter, ale również zwraca ilość wykonywanych podstawień. To
rozwiązanie jest najszybsze z czterech zaprezentowanych. W tym przypadku litery po prostu usuwane, bo
zastępujemy je brakiem znaku.
$a = ($dna =~ tr/Aa//);
$c = ($dna =~ tr/Cc//);
$g = ($dna =~ tr/Gg//);
$t = ($dna =~ tr/Tt//);
Funkcja nie przyjmuje znaków, które nie są literami. Dlatego też ciężko byłoby policzyć znaki, które literami nie są. To
też da się obejść. Najpierw należy usunąć wszystkie zasady, licząc przy okazji ilość usuniętych znaków, a następnie
sprawdzić jaką długość ma pozostały łańcuch.
$zasady = ($dna =~ tr/ACGTacgt//);
$nie_zasady = ((length $dna) - $zasady);
Zadania
1. Napisać program zawierający nieskończoną pętlę. Test warunkowy musi zawsze zwracać prawdę. Aby przerwać
wykonywanie programu można zamknąć okno wiersza poleceń, albo nacisnąć Ctrl + C.
2. Poprosić użytkownika o wprowadzenie dwóch krótkich łańcuchów DNA. Połączyć je ze sobą używając operatora
„.=”. Wydrukować połączone łańcuchy, a następnie wydrukować drugi łańcuch tak, aby otrzymać przedstawioną
poniżej sytuację.
Łańcuchy: AAAA oraz TTTT
AAAATTTT
TTTT
3. Napisać program, który wydrukuje liczby od 1 do 100. Oczywiście proszę użyć mniej niż 100 linii kodu.
4. Napisać program, który wydrukuje sekwencję łańcucha komplementarnego (w orientacji 5’-3’). Proszę nie
korzystać z funkcji s/// lub tr///. Należy użyć funkcji substr i zmieniać zasady pojedynczo.
5. Napisać program do obliczenia procentu hydrofobowych aminokwasów w sekwencji białkowej. Podręcznik do
biochemii może się przydac…
6. Napisać program, który sprawdza, czy dwie sekwencje podane jako argumenty są do siebie komplementarne.
7. Napisać program, który poda procent par GC w danej sekwencji.
8. Napisać program, który zamienia dwie zasady w określonej pozycji w sekwencji DNA. Można spróbować
wykorzystać funkcję substr lub splice.
9. Napisać program, który zapisuje czasowy plik, a następnie usuwa go. Funkcja unlink usuwa plik, na przykład:
unlink „tmpfile”
Program powinien sprawdzić, czy usunięcie pliku powiodło się.