Basis¶

Deze pagina bevat uitvoerbare code.

De vingerafdruk van DNA¶

# Dit notebook staat in solutions/ en de sequenties in problems/assets/dna/,
# vandaar het pad terug. In de opgave zelf is dat gewoon "assets/dna/".
DNA = "../problems/assets/dna/"

database = [
    ["Anouk", 4, 1, 5],
    ["Bram", 5, 2, 3],
    ["Chiara", 3, 5, 1],
]

patterns = ["AGAT", "AATG", "TATC"]

def strip_newline(line):
    """Geeft de regel terug zonder regeleinde aan het eind

    :param line: de regel zoals die uit het bestand komt
    :type line: str
    :rtype: str
    """
    if len(line) > 0 and line[-1] == "\n":
        return line[:-1]

    return line


assert strip_newline("AGAT\n") == "AGAT"
assert strip_newline("AGAT") == "AGAT"
assert strip_newline("") == ""

De controle op de lengte moet eerst. Op een lege string bestaat line[-1] niet en stopt Python met een foutmelding. Omdat and stopt zodra het linkerdeel False is, wordt line[-1] in dat geval nooit opgevraagd.

De laatste regel heeft geen else nodig: als de if toeslaat is de functie via de return al klaar.

def read_sequence(filename):
    """Leest een DNA-sequentie uit een bestand

    :param filename: de naam van het bestand
    :type filename: str
    :rtype: str
    """
    with open(filename) as file:
        for line in file:
            return strip_newline(line)


assert len(read_sequence(DNA + "1.txt")) == 67
assert read_sequence(DNA + "1.txt")[0:4] == "AGCT"
assert read_sequence(DNA + "3.txt")[0:4] == "CCCC"

De return staat in de lus, dus de functie is klaar zodra de eerste regel is gelezen. Dat mag hier, omdat onze sequenties uit één regel bestaan.

Het bestand wordt gesloten zodra het with-blok wordt verlaten, en dat gebeurt ook bij een return van binnenuit. Dat is precies waarvoor die constructie bestaat.

def count_repeats(sequence, pattern, start):
    """Telt hoe vaak het patroon vanaf start direct achter elkaar voorkomt

    :param sequence: de DNA-sequentie
    :type sequence: str
    :param pattern: de STR waarnaar we zoeken
    :type pattern: str
    :param start: de positie waar we beginnen te kijken
    :type start: int
    :rtype: int
    """
    n = 0
    position = start

    while sequence[position:position + len(pattern)] == pattern:
        n = n + 1
        position = position + len(pattern)

    return n


assert count_repeats("AGATAGATCC", "AGAT", 0) == 2
assert count_repeats("AGATAGATCC", "AGAT", 4) == 1
assert count_repeats("AGATAGATCC", "AGAT", 8) == 0
assert count_repeats("CCAGATCC", "AGAT", 0) == 0

Een while-lus, want vooraf weet je niet hoe vaak het patroon zich herhaalt.

Er is hier één ding dat je gratis krijgt: als position voorbij het einde van de sequentie ligt, geeft sequence[position:position + 4] gewoon een kortere of lege string terug in plaats van een foutmelding. Die is dan niet gelijk aan pattern, dus de lus stopt vanzelf. Bij sequence[position] was dat anders geweest.

def longest_match(sequence, pattern):
    """Geeft de langste reeks opeenvolgende herhalingen in de sequentie

    :param sequence: de DNA-sequentie
    :type sequence: str
    :param pattern: de STR waarnaar we zoeken
    :type pattern: str
    :rtype: int
    """
    longest = 0

    for position in range(len(sequence)):
        n = count_repeats(sequence, pattern, position)
        if n > longest:
            longest = n

    return longest


assert longest_match("AGATCCAGATAGATAGAT", "AGAT") == 3
assert longest_match("CCCCCC", "AGAT") == 0
assert longest_match(read_sequence(DNA + "1.txt"), "AGAT") == 4
assert longest_match(read_sequence(DNA + "2.txt"), "TATC") == 3

Elke positie langslopen en de grootste uitkomst onthouden. Dat is niet de snelste aanpak, want je begint ook te tellen op posities waar het patroon helemaal niet staat, maar het is wel de aanpak die je meteen kunt opschrijven en nakijken.

longest begint op 0, en dat is meteen het goede antwoord voor een sequentie waarin het patroon niet voorkomt.

def profile(sequence, patterns):
    """Geeft voor elke STR het aantal herhalingen, in dezelfde volgorde

    :param sequence: de DNA-sequentie
    :type sequence: str
    :param patterns: de STR's waarnaar we kijken
    :type patterns: list
    :rtype: list
    """
    result = []

    for pattern in patterns:
        result = result + [longest_match(sequence, pattern)]

    return result


assert profile(read_sequence(DNA + "1.txt"), patterns) == [4, 1, 5]
assert profile(read_sequence(DNA + "2.txt"), patterns) == [5, 2, 3]
assert profile(read_sequence(DNA + "4.txt"), patterns) == [2, 3, 2]

De lijst wordt opgebouwd met result = result + [...]. De volgorde van patterns bepaalt de volgorde van de uitkomsten, en dat is precies waarom de vergelijking in de volgende stap werkt.

def identify(sequence, patterns, database):
    """Zoekt het profiel van de sequentie op in de database

    :param sequence: de DNA-sequentie
    :type sequence: str
    :param patterns: de STR's waarnaar we kijken
    :type patterns: list
    :param database: de profielen, elk als [naam, aantal, aantal, aantal]
    :type database: list
    :rtype: str
    """
    found = profile(sequence, patterns)

    for row in database:
        if row[1:] == found:
            return row[0]

    return "Geen match"


assert identify(read_sequence(DNA + "1.txt"), patterns, database) == "Anouk"
assert identify(read_sequence(DNA + "2.txt"), patterns, database) == "Bram"
assert identify(read_sequence(DNA + "3.txt"), patterns, database) == "Chiara"
assert identify(read_sequence(DNA + "4.txt"), patterns, database) == "Geen match"

row[1:] is alles behalve de naam, en dat is een lijst van drie getallen. Twee lijsten mag je in één keer met == vergelijken; Python kijkt dan of ze even lang zijn en of alle elementen gelijk zijn.

"Geen match" staat ná de lus. Wordt die regel bereikt, dan is de hele database langsgelopen zonder dat er een return is uitgevoerd, en dus is er niemand gevonden. Zet je hem in de lus, dan geeft de functie al bij de eerste persoon die niet past een antwoord.