Python - wprowadzenie

14 - Zadania 2

Zadanie 1: symulacja zmian częstości alleli w wyniku mutacji

Treść zadania:

Napisz program symulujący zmiany częstości alleli w wyniku mutacji. Zasady i założenia teoretyczne:

  • W populacji dany gen posiada dwa allele, tradycyjnie nazwijmy je: A i a.
  • Allel A mutuje w a z pewnym prawdopodobieństwem, oznaczmy je u.
  • Allel a mutuje w A z innym prawdopodobieństwem (v).
  • W populacji znajduje się początkowo określona liczba jednych i drugich alleli, np. po 100
  • Prawdopodobieństwa mutacji w obie strony i początkowe liczby alleli podaje użytkownik
  • Program ma symulować określoną liczbę pokoleń, w każdym pokoleniu każdy allel A może zmutować w a i odwrotnie.
  • Proporcja alleli w każdym pokoleniu mają być przedstawione liczbowo i ,,symbolicznie’’ (zob. poniżej)
  • Pomiędzy wydrukiem dla kolejnych pokoleń powinna upłynąć sekunda. Znajdź odpowiednią funkcję w module time Standardowej Biblioteki Pythona.
  • Na końcu zapisz wykres przedstawiający zmiany proporcji alleli w czasie.

Przykładowy wydruk podczas działania programu:

Podaj liczbę alleli A: 150
Podaj liczbę alleli a: 50
Podaj p-stwo mutacji A->a (0-1): 0.001
Podaj p-stwo mutacji a->A (0-1): 0.0001
Podaj liczbę pokoleń: 100
Podaj nazwę pliku wykresu [zmiany_alleli.svg]: 
----- Pokolenie: 0, A: 150, a: 50 ------
##################################################

##################################################

##################################################

**************************************************

----- Pokolenie: 1, A: 150, a: 50 ------
##################################################

##################################################

##################################################

**************************************************

----- Pokolenie: 2, A: 149, a: 51 ------
##################################################

##################################################

#################################################*

**************************************************

----- Pokolenie: 3, A: 149, a: 51 ------
##################################################

##################################################

#################################################*

**************************************************

itd...

----- Pokolenie: 999, A: 52, a: 148 ------
##################################################

##************************************************

**************************************************

**************************************************

----- Pokolenie: 1000, A: 52, a: 148 ------
##################################################

##************************************************

**************************************************

**************************************************

Przykładowe rozwiązanie:

Podzielimy kod na funkcje realizujące poszczególne zadania:

  • pobierz_dane() - Pobranie wartości początkowych od użytkownika
  • mutuj() - Symulacja mutacji w pokoleniu
  • wypisz_allele() - Wypisanie alleli w pokoleniu
  • symuluj() - Przeprowadzenie symulacji dla określonej liczby pokoleń (m. in. wywołujemy funkcje mutuj() i wypisz_allele().
  • zrob_wykres() - Wygererowanie, zapisanie i wyświetlenie wykresu.

Następnie wywołamy kolejno pobierz_dane(), symuluj() oraz zrob_wykres(). Funkcje te powinny przekazywać sobie odpowiednie dane.

from random import random
from time import sleep
import matplotlib.pyplot as plt 

def pobierz_dane():
    """Pobranie wartości początkowych od użytkownika"""

    A = int(input('Podaj liczbę alleli A: '))
    a = int(input('Podaj liczbę alleli a: '))
    u = float(input('Podaj p-stwo mutacji A->a (0-1): '))
    v = float(input('Podaj p-stwo mutacji a->A (0-1): '))
    pokolenia = int(input('Podaj liczbę pokoleń: '))
    nazwa_pliku = input('Podaj nazwę pliku wykresu [zmiany_alleli.svg]: ')
    # Jeśli użytkownik nie poda nazwy pliku, to ustawiona jest wartość domyślna
    if nazwa_pliku == '':
        nazwa_pliku = 'zmiany_alleli.svg'
    return A, a, u, v, pokolenia, nazwa_pliku

def mutuj(l_alleli, p_stwo):
    """Symuluje liczbę mutacji dla wystapień danego allelu w populacji."""
    
    zmutowane = 0
    # Dla każdej kopii allelu odbywa się losowanie, czy mutuje
    for i in range(l_alleli):
        # Jeżeli wylosowana została liczba <= p-stwu mutacji, to allel mutuje
        if random() <= p_stwo:
            zmutowane += 1
    return zmutowane

def wypisz_allele(A, a, pokolenie):
    """Wydruk liczby alleli w populacji i w formie symboli. """

    print(f"----- Pokolenie: {pokolenie}, A: {A}, a: {a} ------")
    wszystkie = A + a
    linia = ''
    for i in range(1, wszystkie+1):
        if i <= A:
            linia += '#'
        else:
            linia += '*'
        if i > 0 and i % 50  == 0:
            print(f"{linia}")
            linia = ''


def zrob_wykres(zmiany_A, zmiany_a, nazwa_pliku='zmiany_alleli.svg'):
    """Tworzy wykres zmian liczb alleli."""

    # Linia zielona, ciągła, etykieta: "allel A"   
    plt.plot(zmiany_A, 'g-', label='allel A')
    # Linia niebieska, kropkowana, etykieta: "allel a"  
    plt.plot(zmiany_a, 'r-', label='allel b')
    plt.xlabel('pokolenia',fontsize=15)
    plt.ylabel('liczba alleli',fontsize=15)
    # legenda loc=0 oznacza optymalną lokalizację
    plt.legend(loc=0)
    plt.savefig(nazwa_pliku)

def symuluj(A, a, u, v, pokolenia):
    """Symulacja zmian liczb alleli."""

    zmiany_A = []
    zmiany_a = []
    for i in range(pokolenia+1):
        wypisz_allele(A, a, i)
        zmiana_A = mutuj(A, u)
        zmiana_a = mutuj(a, v)
        A = A - zmiana_A + zmiana_a
        a = a - zmiana_a + zmiana_A
        zmiany_A.append(A)
        zmiany_a.append(a)
        # Wstrzymanie działania na 1 sek.
        sleep(1)
    return zmiany_A, zmiany_a

# Pobranie danych
A, a, u, v, pokolenia, nazwa_pliku = pobierz_dane()
# Symulacja
zmiany_A, zmiany_a = symuluj(A, a, u, v, pokolenia)
# Tworzenie wykresu
zrob_wykres(zmiany_A, zmiany_a, nazwa_pliku)

Zadanie 2: Symulacja niezależnej ewolucji dwu, początkowo identycznych sekwencji - zmiany dystansu ewolucyjnego między nimi.

Od momentu specjacji, czyli kiedy z jednego gatunku powstaną dwa, ewolucja sekwencji DNA biegnie w nich niezależnie. Zatem z czasem, w skutek kumulacji mutacji, sekwencje DNA, nawet jeśli początkowo są identyczne, różnią się od siebie coraz bardziej. Różnice ewolucyjne między sekwencjami można określić obliczając dystans ewolucyjny/odległość ewolucyjną między nimi. Jest wiele metod liczenia takiego dystansu, najprostsze opierają się na obliczeniu liczby różnic między sekwencjami w stosunku do długości sekwencji. Więcej na ten temat można poczytać np. na stronie: https://ggoralski.gitbook.io/filogenetyka/teoria/04-modele_ewolucji_molekularnej.

Program, który należy napisać w tym zadaniu, powinien pokazać w jaki sposób zmienia się dystans ewolucyjny między dwoma, początkowo identycznymi sekwencjami w czasie. Podczas symulacji sekwencje mutują, przy czym uwzględniamy wyłącznie substytucje, czyli zmianę jednej zasady w inną. Dla każdego pokolenia obliczany jest dystans ewolucyjny między sekwencjami a na końcu tworzony jest wykres pokazujący zmiany dystansu w trakcie symulacji.

Założenia i zasada działania:

  • Użytkownik podaje:
    • Początkową sekwencję DNA, lub jej długość. W drugim przypadku tworzony jest losowy ciąg zasad o podanej długości.
    • Prawdopodobieństwo mutacji/zasadę w pokoleniu
    • Liczbę pokoleń
    • Nazwę pliku wykresu - w przypadku niepodania, używana jest nazwa domyślna.
  • Symulacja trwa przez podaną przez użytkownika liczbę ,,pokoleń’’. W każdym pokoleniu:
    • Obliczany jest dystans ewolucyjny między sekwencjami - w tym przypadku jest to liczba różnic między odpowiadającymi sobie (homologicznymi) miejscami w sekwencjach podzielona przez całkowitą liczbę zasad sekwencji. Jest on zapamiętywany do stworzenia wykresu na końcu działania programu.
    • W terminalu drukowane są:
      • Pokolenie
      • Dystans
      • Obie sekwencje
    • Losowane jest zaistnienie mutacji dla każdej zasady w obu sekwencjach. W przypadku wylosowania mutacji, zasada zmienia się na losową inną zasadę, np. A na C lub G lub T.
  • Po zakończeniu symulacji tworzony jest wykres wartości dystansu ewolucyjnego miedzy sekwencjami w czasie symulacji, następnie jest on spisywany w pliku i wyświetlany.
  • Poszczególne zadania powinny być umieszczone w oddzielnych funkcjach, można też utworzyć osobny moduł(y).

Przykładowe rozwiązanie:

Dzielimy kod na funkcje realizujące poszczególne zadania:

  • pobierz_ustawienia() - Pobranie ustawień symulacji od użytkownika.
  • oblicz_dystans() - Obliczenie dystansu między dwiema sekwencjami.
  • mutuj() - Generowanie mutacji dla sekwencji w danym pokoleniu.
  • symuluj() - Symulacja przez określoną liczbę pokoleń. Wywołuje oblicz_dystans() i mutuj().
  • generuj_wykres() - Wygenerowanie, zapisanie i wyświetlenie wykresu.

Następnie wywołujemy kolejno: pobierz_ustawienia(), symuluj() i generuj_wykres().

import random as rnd 
import matplotlib.pyplot as plt 

def pobierz_ustawienia():
    """ Pobiera wartości początkowe symulacji od użytkownika."""

    sekwencja = input('Podaj sekwencję wyjściową, jeśli chcesz losować, wciśnij Enter: ')
    # Jeśli uzytkownik nie poda sekwencji, generowana jest 
    # sekwencja losowa o podanej długości
    if len(sekwencja) == 0:
        dlugosc = int(input('Podaj długość sekwencji: '))
        zasady = 'ACGT'
        sekwencja = rnd.choices(zasady, k=dlugosc)

    pstwo_mutacji = float(input('Podaj p-stwo mutacji/zasadę w pokoleniu: '))
    l_pokolen = int(input('Podaj liczbę pokoleń: '))
    nazwa_pliku = input('Podaj nazwę pliku wykresu [dystanse.svg]: ')
    if nazwa_pliku == '':
        nazwa_pliku = 'dystanse.svg'
    return sekwencja, pstwo_mutacji, l_pokolen, nazwa_pliku

def oblicz_dystans(sekw_1, sekw_2):
    """Oblicza dystans między dwoma sekwencjami.
    
       Dystans jest liczony jako stosunek liczby róznic do długości sekwencji.
    """
    roznice = 0
    # Porównanie obu sekwencji, znak po znaku.
    for l_1, l_2 in zip(sekw_1, sekw_2):
        if l_1 != l_2:
            roznice += 1
    dystans = roznice / len(sekw_1)
    return dystans

def mutuj(sekwencja, pstwo_mutacji):
    """Symuluje mutacje danej sekwencji w pokoleniu."""

    # Słownik z możliwymi mutacjami dla zasad:
    mutacje = {
        'A': 'CGT',
        'T': 'CGA',
        'C': 'AGT',
        'G': 'ACT'
    }
    # Dla każdej zasady w sekwencji jest możliwa mutacja
    for i in range(len(sekwencja)):
            # Jeżeli wylosowana liczba jest <= ps-stwu mutacji
            # to zachodzi mutacja
            if rnd.random() <= pstwo_mutacji:
                # Zasada w danym miejscu zmienia się na losową inną zasadę
                sekwencja[i] = rnd.choice(mutacje[sekwencja[i]])
    return sekwencja

def symuluj(sekwencja, pstwo_mutacji, l_pokolen):
    """Symuluje niezależną ewolucję dwu sekwencji, pochodzących od jednej."""

    # Tworzymy dwie listy reprezentujące dwie takie same sekwencje.
    sekw_1 = list(sekwencja)
    sekw_2 = list(sekwencja)
    # Lista zachowująca kolejne dystanse między sekwencjami
    dystanse = []
    # Symulacja przez l_pokolen
    for pokolenie in range(l_pokolen+1):
        dystans = oblicz_dystans(sekw_1, sekw_2)
        dystanse.append(dystans)
        print(f"Pokolenie: {pokolenie}, dystans: {dystans}")
        # Drukowanie bieżących sekwencji
        print(''.join(sekw_1))
        print(''.join(sekw_2))
        sekw_1 = mutuj(sekw_1, pstwo_mutacji)
        sekw_2 = mutuj(sekw_2, pstwo_mutacji)
    return(dystanse)

def generuj_wykres(dystanse, nazwa_pliku):
    """Generuje wykres zmian dystansu między sekwencjami."""

    # Tworzenie wykresu    
    plt.plot(dystanse)
    # etykiety osi x i y
    plt.xlabel('pokolenie',fontsize=15)
    plt.ylabel('dystans',fontsize=15)
    # Zapisanie pliku w formacie svg
    plt.savefig(nazwa_pliku)
    # Wyświetlenie wykresu
    plt.show()

# Pobranie ustawień początkowych
sekwencja, pstwo_mutacji, l_pokolen, nazwa_pliku = pobierz_ustawienia()
# Przebieg symulacji, pobranie zmian dystansu miedzy sekwencjami
dystanse = symuluj(sekwencja, pstwo_mutacji, l_pokolen)
# Generowanie wykresu
generuj_wykres(dystanse, nazwa_pliku)

Podaj sekwencję wyjściową, jeśli chcesz losować, wciśnij Enter:  
Podaj długość sekwencji:  80
Podaj p-stwo mutacji/zasadę w pokoleniu:  0.05
Podaj liczbę pokoleń:  200
Podaj nazwę pliku wykresu [dystanse.svg]:  

Pokolenie: 0, dystans: 0.0
ACAGGCCAGCAACATGCTTTTGGCAGTCCTGAGTGGAACGCCCACGGGAATCTTGTCACAGATGTAGGTCATTATATACA
ACAGGCCAGCAACATGCTTTTGGCAGTCCTGAGTGGAACGCCCACGGGAATCTTGTCACAGATGTAGGTCATTATATACA
Pokolenie: 1, dystans: 0.125
ACAGGCCAGCAACATGCTTTTTGCAGTCCTGAGAGGAACGCCCACGGGAATCTTGTCACAGATGTGTGTCATTATATACA
ACGGGCCAGCAACATGCTTTTGGCAGTCATGAGTGGAACGCCAACGCGAATCTTGTCACAGATCGAGGTCATTATATACA

...

Pokolenie: 199, dystans: 0.775
GTGAGGTCTCACGCAGGAATGTGGGTAATACTGGCCTGTATTTAAAAGGCAACACGTTATACGGATTTGAGTTTCATTGC
CTCATTTAGGGCGTGGCGCTACCACGCTAATACATAACTCCCAATTCGATGGGTTCTATTCTAATGGACAATGACTGTAC
Pokolenie: 200, dystans: 0.75
GTGAGGTCTCACGCAGGAATGTGGGTAATGCTGGCCTGTACTTAAAAGGCAACACGTTATACGGATTTGAGTTTCATTGC
CTCATTTAGGGCGTGGCGCTACCACGATAATACAGAACTCCCAATTCGATGGGTTCTATTCTAATGGAGAATGACTGTAC

Czyszczenie ekranu

Powyższe programy wyświetlają wyniki dla kolejnych pokoleń w formie kolejnych linii znaków, co nie sprzyja estetyce i łatwości obserwacji zmian. Można to nieco zmienić, czyszcząc ekran terminala między pokoleniami co pozwala uzyskać efekt zmieniających się danych albo nawet ,,animacji".

Utworzymy w tym celu funkcję ,,czyszczącą" ekran. Pewien problem tkwi w tym, że nieco inaczej wygląda to w systemie Windows a inaczej w systemach takich jak GNU Linux, macOS itd.. Na szczęście moduł os Biblioteki Standardowej Pythona pozwala wykryć rodzaj systemu operacyjnego i wywołać odpowiednie polecenie właściwe dla systemu. Wywołując os.name uzyskujemy nt (dla Windows), posix (Linux, macOS itd) lub java. Funkcja os.system() pozwala przekazać polecenie do powłoki systemu (tak, jakbyśmy wpisali je w terminalu).

import os

def wyczysc(): 
  
    # Jeśli windows, os.name zwraca 'nt'
    if os.name == 'nt': 
    	# Czyszczenie ekranu w systemie Windows
        os.system('cls') 
  
    # W innym przypadku otrzymamy 'posix'
    else: 
        # Czyszczenie ekranu w Linukxie, macOS itd
        os.system('clear')

Funkcję wyczyść() zwykle najlepiej wywołać w kodzie bezpośrednio przed wydrukowaniem na ekranie bieżących wartości. Można też dodatkowo użyć funkcji sleep() z modułu time, aby zmiany nie były zbyt szybkie.

Zmodyfikuj powyższe programy używając funkcji wyczysc().

Last updated on 15 Dec 2020
Published on 15 Dec 2020