diff --git a/BINF_TP5.ipynb b/BINF_TP5.ipynb
index 9a6c7bf..40e607c 100644
--- a/BINF_TP5.ipynb
+++ b/BINF_TP5.ipynb
@@ -1,269 +1,369 @@
{
- "nbformat": 4,
- "nbformat_minor": 0,
- "metadata": {
- "colab": {
- "provenance": [],
- "authorship_tag": "ABX9TyO0SUJI6kaczFcOh8NoKcqb",
- "include_colab_link": true
- },
- "kernelspec": {
- "name": "python3",
- "display_name": "Python 3"
- },
- "language_info": {
- "name": "python"
- }
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "colab_type": "text",
+ "id": "view-in-github"
+ },
+ "source": [
+ "
"
+ ]
},
- "cells": [
- {
- "cell_type": "markdown",
- "metadata": {
- "id": "view-in-github",
- "colab_type": "text"
- },
- "source": [
- "
"
- ]
- },
- {
- "cell_type": "markdown",
- "source": [
- "# *TP5 Reconstruction de génome par alignement sur une référence - Algorithme bowtie et transformée de Burrows-Wheeler*\n",
- "\n",
- "On cherche à reconstruire le génome d'un organisme dont l'espèce a déjà été séquencé. On va utiliser le génome existant pour guider notre reconstruction en cherchant la position de chaque fragment séquencé dans le génome de référence.\n",
- "\n",
- "On appelle read un fragment séquencé. La taille d'un read dépend de la mêthode de séquençage. La méthode de Sanger produit des reads d'environ 1000 nucléotides alors que les méthodes de séquençage nouvelle génération vont avoir des reads plus courts d'environ 100 nucléotides.\n",
- "\n",
- "Les génomes ayant des tailles de l'ordre de 100 000 paires de bases et plus, leur séquençage requiert l'acquisition de milliers voirs millions de reads dont il faut trouver la position dans le génome de référence. Pour que cette méthode soit viable, il faut être capable d'effectuer cette recherche de manière efficace.\n",
- "\n",
- "La méthode bowtie propose une solution à ce problème basé sur le calcul de la tranformée de Burrows-Wheeler du génome de référence. La transformée de Burrow-Wheeler est une structure basée sur l'ordonancement des préfixes / suffixes d'une séquence permettant une recherche efficace de sous-séquence. Pour pouvoir faire le lien entre transformée et séquence d’origine, on va aussi calculer un index faisant la correspondance entre un caractère de la séquence transformée et la position d’apparition du préfixe correspondant.\n",
- "\n",
- ""
- ],
- "metadata": {
- "id": "BoOWxsamkRHR"
- }
- },
- {
- "cell_type": "markdown",
- "source": [
- "## Exercice 1 : Transformée de Burrows-Wheeler\n",
- "\n",
- "La transformée de Burrows-Wheeler (BWT) s’obtient en triant par ordre lexicographique toutes les permutations circulaires possibles d’une séquence, et en conservant le dernier caractère de chacune.\n",
- "\n",
- "Ecrire la fonction BWT(s) -> str qui retourne la transformée BW de la séquence s en effectuant les étapes suivantes :\n",
- "\n",
- "\n",
- "1. Ajouter un caractère $ à la fin de s (celui-ci servira pour la transformée inverse).\n",
- "2. Générer la liste perm toutes les permutations circulaires de s (on fait passer la dernière lettre au début de s).\n",
- "3. Générer un tableau d’indice qui contient l’ordre de chaque permutation.\n",
- "4. Trier le tableau d’indices et perm par ordre lexicographique de perm.\n",
- "5. Retourner la séquence contenant le dernier caractère de chaque permutation dans l’ordre lexicographique ainsi que le tableau d’indices.\n"
- ],
- "metadata": {
- "id": "Wc6rm1agmoCb"
- }
- },
- {
- "cell_type": "code",
- "source": [
- "print(\"votre fonction ici !!\")"
- ],
- "metadata": {
- "id": "D5zJsnmQnajj"
- },
- "execution_count": null,
- "outputs": []
- },
- {
- "cell_type": "markdown",
- "source": [
- "### Exemple\n",
- "\n",
- "Transformée de chien :\n",
- "\n",
- "* Ajouter un caractère \\\\$ à la fin : chien$\n",
- "\n",
- "* Générer le tableau perm des permutations cirulaires :\n",
- "\n",
- "$$\n",
- "\\begin{array}{l}\n",
- "chien$\\\\\n",
- "$chien\\\\\n",
- "n$chie\\\\\n",
- "en$chi\\\\\n",
- "ien$ch\\\\\n",
- "hien$c\n",
- "\\end{array}\n",
- "$$\n",
- "* Trier ce tableau par ordre lexicographique :\n",
- "$$\n",
- "\\begin{array}{l}\n",
- "$chien\\\\\n",
- "chien$\\\\\n",
- "en$chi\\\\\n",
- "hien$c\\\\\n",
- "ien$ch\\\\\n",
- "n$chie\n",
- "\\end{array}\n",
- "$$\n",
- "* La transformée est la dernière colonne de ce tableau (le dernier caractère pour chaque ligne) : n$iche\n"
- ],
- "metadata": {
- "id": "bG6dMzI1ndK-"
- }
- },
- {
- "cell_type": "markdown",
- "source": [
- "## Exercice 2 : Transformée inverse\n",
- "\n",
- "Pour vérifier que la transformée est correcte, on va implémenter la transformée inverse. Pour ce faire, on va recréer l’ensemble des permutations de manière itérative en faisant un nombre d’itérations égal à la taille de la transformée. On va commencer avec un ensemble de séquences vides qu’on appelle table de même taille que la transformée qu’on va modifier à chaque itération :\n",
- "\n",
- "* Pour chaque caractère en position i de la transformee, l’ajouter au début de la ième séquence de table.\n",
- "* Trier table par ordre lexicographique.\n",
- "\n",
- "Le résultat est la séquence de table qui se termine par le caractère $. Notez que c’est pour ça qu’on a besoin d’ajouter ce caractère au moment de la transformée.\n",
- "\n",
- "Ecrivez la fonction iBWT(t) -> str qui retourne la transformée inverse de la transformée BW t.\n"
- ],
- "metadata": {
- "id": "JgSWho3GorwJ"
- }
- },
- {
- "cell_type": "code",
- "source": [
- "print(\"Votre code ici !\")"
- ],
- "metadata": {
- "id": "oH-bAVbgpOUe"
- },
- "execution_count": null,
- "outputs": []
- },
- {
- "cell_type": "markdown",
- "source": [
- "### Exemple :\n",
- "Transformée inverse de n$iche :\n",
- "\n",
- "$$\n",
- "\\begin{array}{lllllllllll}\n",
- "ADD1 & SORT1 & ADD2 & SORT2 & ADD3 & SORT3 & ADD4 & SORT4 & ADD5 & SORT5 & ADD6\\\\\n",
- "n & $ & n$ & $c & n$c & $ch & n$ch & $chi & n$chi & $chie & n$chie\\\\\n",
- "$ & c & $c & ch & $ch & chi & $chi & chie & $chie & chien & $chien\\\\\n",
- "i & e & ie & en & ien & en$ & ien$ & en$c & ien$c & en$ch & ien$ch\\\\\n",
- "c & h & ch & hi & chi & hie & chie & hien & chien & hien$ & chien$\\\\\n",
- "h & i & hi & ie & hie & ien & hien & ien$ & hien & ien$c & ien$c\\\\\n",
- "e & n & en & n$ & en$ & n$c & en$c & n$ch & en$ch & n$chi & en$chi\\\\\n",
- "\\end{array}\n",
- "$$\n"
- ],
- "metadata": {
- "id": "MMRvEACHpVDs"
- }
- },
- {
- "cell_type": "markdown",
- "source": [
- "### Exercice 3 : Recherche d’un mot dans la transformée\n",
- "\n",
- "Maintenant que l’on dispose de la transformée, on va développer un algorithme permettant de trouver efficacement si un mot se trouve dedans. Pour cela on va implémenter l’algorithme suivant qui provient de l’outil Bowtie base sur un principe similaire a la construction de table de l’exercice précédent. L’algorithme va reconstruire itérativement l’ensemble des lignes de la table contenant des fragments de plus en plus longs du mot recherche, jusqu’à obtenir l’ensemble des lignes contenant le mot entier. Vu que la table est toujours triée lexicographiquement, on sait que ces séquences vont se trouver cote a cote. L’algorithme est un peu plus complexe que la transformée inverse de l’exercice précédent car il ne reconstruit pas la table entière mais se base sur les propriétés de celle-ci pour trouver les indices :\n",
- "\n",
- "\n",
- "\n",
- "Ou P est la séquence recherchée de taille p, C est un tableau qui donne pour chaque caractère le nombre d’occurrences dans la transformée de caractères qui lui sont lexicographiquement inferieurs et Occ(c,k) est une fonction qui retourne le nombre d’occurrences du caractère c dans la transformée jusqu’à la position k (non-inclue). Les indices sp et ep indiquent le début et la fin respectivement des lignes contenant un fragment de P dans la table. Si P est présente dans la transformée, alors elle sera présente dans toutes les séquences entre les indices sp et ep. Le tableau C peut être précalculé. Attention dans l’algorithme ci-dessus les indices commencent à 1.\n",
- "\n",
- "Ecrivez la fonction exactmatch(p, bwt) -> (int, int) qui retourne la paire d'indices (sp, ep) (avec sp < ep) de toutes les séquences dans la transformée bwt qui contiennent la sous-séquence s."
- ],
- "metadata": {
- "id": "ZveHazSWsOJ_"
- }
- },
- {
- "cell_type": "code",
- "source": [
- "print(\"Votre code ici !\")"
- ],
- "metadata": {
- "id": "C63yJjfjtUSQ"
- },
- "execution_count": null,
- "outputs": []
- },
- {
- "cell_type": "markdown",
- "source": [
- "## Exercice 4: retrouver les indices dans la séquence d'origine\n",
- "\n",
- "A ce stade, il reste une feinte !! En effet, les indices que l’on a récupérés sont des indices dans la transformée et non dans la séquence d’origine. Or ce qui nous intéresse c’est d’obtenir toutes les positions ou apparaissent la séquence recherchée dans la séquence d’origine. Pour ce faire, on va se baser sur le tableau d’indices que l’on a calculé dans la transformée. Ce tableau fait un lien régulier entre une position dans la transformée et position dans la séquence de base. Pour chaque position sp≤k≤se on va chercher si k est dans le tableau d’indice, si oui alors on le retourne, sinon on va dérouler la transformée inverse jusqu’à tombe sur un indice présent dans le tableau d’indice et on retourne la position correspondante plus le nombre d’étapes de transformées inverses effectués. L’algorithme est le suivant :\n",
- "\n",
- "\n",
- "\n",
- "ou k est l’indice recherché, idxs le tableau d’indices et stepleft est la fonction suivante :\n",
- "\n",
- "\n",
- "\n",
- "Notez que le modulo apparait car on peut se retrouver à dépasser la taille de la transformée initiale dans la recherche de l’indice le plus proche.\n",
- "\n",
- "Ecrivez la fonction index_from_match(k, idxs) -> int qui retourne l'indice dans la séquence d'origine à partir de k l'indice dans la transformée et idxs le tableau des correspondances des positions."
- ],
- "metadata": {
- "id": "uRHVKu1s9n1l"
- }
- },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "BoOWxsamkRHR"
+ },
+ "source": [
+ "# *TP5 Reconstruction de génome par alignement sur une référence - Algorithme bowtie et transformée de Burrows-Wheeler*\n",
+ "\n",
+ "On cherche à reconstruire le génome d'un organisme dont l'espèce a déjà été séquencé. On va utiliser le génome existant pour guider notre reconstruction en cherchant la position de chaque fragment séquencé dans le génome de référence.\n",
+ "\n",
+ "On appelle read un fragment séquencé. La taille d'un read dépend de la mêthode de séquençage. La méthode de Sanger produit des reads d'environ 1000 nucléotides alors que les méthodes de séquençage nouvelle génération vont avoir des reads plus courts d'environ 100 nucléotides.\n",
+ "\n",
+ "Les génomes ayant des tailles de l'ordre de 100 000 paires de bases et plus, leur séquençage requiert l'acquisition de milliers voirs millions de reads dont il faut trouver la position dans le génome de référence. Pour que cette méthode soit viable, il faut être capable d'effectuer cette recherche de manière efficace.\n",
+ "\n",
+ "La méthode bowtie propose une solution à ce problème basé sur le calcul de la tranformée de Burrows-Wheeler du génome de référence. La transformée de Burrow-Wheeler est une structure basée sur l'ordonancement des préfixes / suffixes d'une séquence permettant une recherche efficace de sous-séquence. Pour pouvoir faire le lien entre transformée et séquence d’origine, on va aussi calculer un index faisant la correspondance entre un caractère de la séquence transformée et la position d’apparition du préfixe correspondant.\n",
+ "\n",
+ ""
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "Wc6rm1agmoCb"
+ },
+ "source": [
+ "## Exercice 1 : Transformée de Burrows-Wheeler\n",
+ "\n",
+ "La transformée de Burrows-Wheeler (BWT) s’obtient en triant par ordre lexicographique toutes les permutations circulaires possibles d’une séquence, et en conservant le dernier caractère de chacune.\n",
+ "\n",
+ "Ecrire la fonction BWT(s) -> str qui retourne la transformée BW de la séquence s en effectuant les étapes suivantes :\n",
+ "\n",
+ "\n",
+ "1. Ajouter un caractère $ à la fin de s (celui-ci servira pour la transformée inverse).\n",
+ "2. Générer la liste perm toutes les permutations circulaires de s (on fait passer la dernière lettre au début de s).\n",
+ "3. Générer un tableau d’indice qui contient l’ordre de chaque permutation.\n",
+ "4. Trier le tableau d’indices et perm par ordre lexicographique de perm.\n",
+ "5. Retourner la séquence contenant le dernier caractère de chaque permutation dans l’ordre lexicographique ainsi que le tableau d’indices.\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 24,
+ "metadata": {
+ "id": "D5zJsnmQnajj"
+ },
+ "outputs": [
{
- "cell_type": "code",
- "source": [
- "print(\"Votre code ici !\")"
- ],
- "metadata": {
- "id": "3RzUM2FJ91u_"
- },
- "execution_count": null,
- "outputs": []
- },
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "['$chien', 'chien$', 'en$chi', 'hien$c', 'ien$ch', 'n$chie']\n",
+ "n$iche [1 0 3 5 4 2]\n"
+ ]
+ }
+ ],
+ "source": [
+ "import numpy as np\n",
+ "\n",
+ "def BWT(s):\n",
+ " curr = s + \"$\"\n",
+ " perm = [curr]\n",
+ " for i in range(len(s)):\n",
+ " curr = curr[-1] + curr[:-1]\n",
+ " perm.append(curr)\n",
+ "\n",
+ " inds = np.argsort(perm)\n",
+ " perm.sort()\n",
+ " print(perm)\n",
+ " seq = \"\"\n",
+ " for w in perm:\n",
+ " seq += w[-1]\n",
+ " return seq, inds\n",
+ "t, inds = BWT(\"chien\")\n",
+ "print(t, inds)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "bG6dMzI1ndK-"
+ },
+ "source": [
+ "### Exemple\n",
+ "\n",
+ "Transformée de chien :\n",
+ "\n",
+ "* Ajouter un caractère \\\\$ à la fin : chien$\n",
+ "\n",
+ "* Générer le tableau perm des permutations cirulaires :\n",
+ "\n",
+ "$$\n",
+ "\\begin{array}{l}\n",
+ "chien$\\\\\n",
+ "$chien\\\\\n",
+ "n$chie\\\\\n",
+ "en$chi\\\\\n",
+ "ien$ch\\\\\n",
+ "hien$c\n",
+ "\\end{array}\n",
+ "$$\n",
+ "* Trier ce tableau par ordre lexicographique :\n",
+ "$$\n",
+ "\\begin{array}{l}\n",
+ "$chien\\\\\n",
+ "chien$\\\\\n",
+ "en$chi\\\\\n",
+ "hien$c\\\\\n",
+ "ien$ch\\\\\n",
+ "n$chie\n",
+ "\\end{array}\n",
+ "$$\n",
+ "* La transformée est la dernière colonne de ce tableau (le dernier caractère pour chaque ligne) : n$iche\n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "JgSWho3GorwJ"
+ },
+ "source": [
+ "## Exercice 2 : Transformée inverse\n",
+ "\n",
+ "Pour vérifier que la transformée est correcte, on va implémenter la transformée inverse. Pour ce faire, on va recréer l’ensemble des permutations de manière itérative en faisant un nombre d’itérations égal à la taille de la transformée. On va commencer avec un ensemble de séquences vides qu’on appelle table de même taille que la transformée qu’on va modifier à chaque itération :\n",
+ "\n",
+ "* Pour chaque caractère en position i de la transformee, l’ajouter au début de la ième séquence de table.\n",
+ "* Trier table par ordre lexicographique.\n",
+ "\n",
+ "Le résultat est la séquence de table qui se termine par le caractère $. Notez que c’est pour ça qu’on a besoin d’ajouter ce caractère au moment de la transformée.\n",
+ "\n",
+ "Ecrivez la fonction iBWT(t) -> str qui retourne la transformée inverse de la transformée BW t.\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 28,
+ "metadata": {
+ "id": "oH-bAVbgpOUe"
+ },
+ "outputs": [
{
- "cell_type": "markdown",
- "source": [
- "## Exercice 5 : Transformée progressive\n",
- "\n",
- "Le calcul de la transformée BW a un énorme coût en mémoire (il faut déterminer l’ordre des préfixes pour chaque position) ce qui la rend inutilisable pour des génomes. Une version légèrement modifiée permet d’échanger de la consommation mémoire contre du temps de calcul en séparant l’ensemble des suffixes à trier en k sous-ensembles, générer la transformée pour chaque sous-ensemble et les concaténer. Pour générer les sous-ensembles, on va simplement tirer k positions aléatoires le long de la séquence.\n",
- "\n",
- "Ecrivez la fonction BWT_prog(s, k) -> str qui retourne la transformée de la séquence s en utilisant k graînes.\n"
- ],
- "metadata": {
- "id": "PtBlbmeNtaA2"
- }
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "['$', 'c', 'e', 'h', 'i', 'n']\n",
+ "['$c', 'ch', 'en', 'hi', 'ie', 'n$']\n",
+ "['$ch', 'chi', 'en$', 'hie', 'ien', 'n$c']\n",
+ "['$chi', 'chie', 'en$c', 'hien', 'ien$', 'n$ch']\n",
+ "['$chie', 'chien', 'en$ch', 'hien$', 'ien$c', 'n$chi']\n",
+ "['$chien', 'chien$', 'en$chi', 'hien$c', 'ien$ch', 'n$chie']\n"
+ ]
},
{
- "cell_type": "code",
- "source": [
- "print(\"Votre code ici !\")"
- ],
- "metadata": {
- "id": "aoWF2thh8zEn"
- },
- "execution_count": null,
- "outputs": []
- },
+ "data": {
+ "text/plain": [
+ "'chien$'"
+ ]
+ },
+ "execution_count": 28,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "def iBWT(t):\n",
+ " table = [\"\" for i in range(len(t))]\n",
+ " for _ in range(len(t)):\n",
+ " for i in range(len(t)):\n",
+ " table[i] = t[i] + table[i]\n",
+ " table.sort()\n",
+ " print(table)\n",
+ " for w in table:\n",
+ " if w[-1] == '$':\n",
+ " return w\n",
+ "iBWT(t)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "MMRvEACHpVDs"
+ },
+ "source": [
+ "### Exemple :\n",
+ "Transformée inverse de n$iche :\n",
+ "\n",
+ "$$\n",
+ "\\begin{array}{lllllllllll}\n",
+ "ADD1 & SORT1 & ADD2 & SORT2 & ADD3 & SORT3 & ADD4 & SORT4 & ADD5 & SORT5 & ADD6\\\\\n",
+ "n & $ & n$ & $c & n$c & $ch & n$ch & $chi & n$chi & $chie & n$chie\\\\\n",
+ "$ & c & $c & ch & $ch & chi & $chi & chie & $chie & chien & $chien\\\\\n",
+ "i & e & ie & en & ien & en$ & ien$ & en$c & ien$c & en$ch & ien$ch\\\\\n",
+ "c & h & ch & hi & chi & hie & chie & hien & chien & hien$ & chien$\\\\\n",
+ "h & i & hi & ie & hie & ien & hien & ien$ & hien & ien$c & ien$c\\\\\n",
+ "e & n & en & n$ & en$ & n$c & en$c & n$ch & en$ch & n$chi & en$chi\\\\\n",
+ "\\end{array}\n",
+ "$$\n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "ZveHazSWsOJ_"
+ },
+ "source": [
+ "### Exercice 3 : Recherche d’un mot dans la transformée\n",
+ "\n",
+ "Maintenant que l’on dispose de la transformée, on va développer un algorithme permettant de trouver efficacement si un mot se trouve dedans. Pour cela on va implémenter l’algorithme suivant qui provient de l’outil Bowtie base sur un principe similaire a la construction de table de l’exercice précédent. L’algorithme va reconstruire itérativement l’ensemble des lignes de la table contenant des fragments de plus en plus longs du mot recherche, jusqu’à obtenir l’ensemble des lignes contenant le mot entier. Vu que la table est toujours triée lexicographiquement, on sait que ces séquences vont se trouver cote a cote. L’algorithme est un peu plus complexe que la transformée inverse de l’exercice précédent car il ne reconstruit pas la table entière mais se base sur les propriétés de celle-ci pour trouver les indices :\n",
+ "\n",
+ "\n",
+ "\n",
+ "Ou P est la séquence recherchée de taille p, C est un tableau qui donne pour chaque caractère le nombre d’occurrences dans la transformée de caractères qui lui sont lexicographiquement inferieurs et Occ(c,k) est une fonction qui retourne le nombre d’occurrences du caractère c dans la transformée jusqu’à la position k (non-inclue). Les indices sp et ep indiquent le début et la fin respectivement des lignes contenant un fragment de P dans la table. Si P est présente dans la transformée, alors elle sera présente dans toutes les séquences entre les indices sp et ep. Le tableau C peut être précalculé. Attention dans l’algorithme ci-dessus les indices commencent à 1.\n",
+ "\n",
+ "Ecrivez la fonction exactmatch(p, bwt) -> (int, int) qui retourne la paire d'indices (sp, ep) (avec sp < ep) de toutes les séquences dans la transformée bwt qui contiennent la sous-séquence s."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 35,
+ "metadata": {
+ "id": "C63yJjfjtUSQ"
+ },
+ "outputs": [
{
- "cell_type": "markdown",
- "source": [
- "## Exercice 6: reconstruction de génome\n",
- "\n",
- "Téléchargez le génome de la mitochondrie humaine: https://www.dropbox.com/scl/fi/ux5m9r80ejm4lz7w3qv0c/sequence_mito.fasta?rlkey=qoaz38vwmg3fgkymx09c7m59j&dl=0\n",
- "Téléchargez l'ensemble des reads ici: https://www.dropbox.com/scl/fi/uajxylv3yfo1is876x1e9/reads_sequence_mito_l69_seed43_N2000_mu0_norc.fasta?rlkey=h16uoh8xjostqbbmy1tc1k86d&dl=0\n",
- "\n",
- "1. Calculez la transformée BW du génome\n",
- "2. Alignez chaque read sur le génome\n",
- "3. Calculez la séquence concessus"
- ],
- "metadata": {
- "id": "RsY4KN-x83cM"
- }
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "{'n': 5, '$': 0, 'i': 4, 'c': 1, 'h': 3, 'e': 2}\n"
+ ]
}
- ]
-}
\ No newline at end of file
+ ],
+ "source": [
+ "def Occ(c, k, bwt):\n",
+ " res = 0\n",
+ " for i in range(k):\n",
+ " if bwt[i] == c:\n",
+ " res += 1\n",
+ " return res\n",
+ "\n",
+ "def C(bwt):\n",
+ " res = {}\n",
+ " for c in bwt:\n",
+ " nbinfs = 0\n",
+ " for car in bwt:\n",
+ " if car < c:\n",
+ " nbinfs += 1\n",
+ " res[c] = nbinfs\n",
+ " return res\n",
+ "#test de C\n",
+ "#print(C(\"n$iche\"))\n",
+ "\n",
+ "def exactmatch(p, bwt):\n",
+ " C = C(p, bwt)\n",
+ " c = p[-1]\n",
+ " sp = C[c] + 1\n",
+ " ep = C[c+1] + 1"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "uRHVKu1s9n1l"
+ },
+ "source": [
+ "## Exercice 4: retrouver les indices dans la séquence d'origine\n",
+ "\n",
+ "A ce stade, il reste une feinte !! En effet, les indices que l’on a récupérés sont des indices dans la transformée et non dans la séquence d’origine. Or ce qui nous intéresse c’est d’obtenir toutes les positions ou apparaissent la séquence recherchée dans la séquence d’origine. Pour ce faire, on va se baser sur le tableau d’indices que l’on a calculé dans la transformée. Ce tableau fait un lien régulier entre une position dans la transformée et position dans la séquence de base. Pour chaque position sp≤k≤se on va chercher si k est dans le tableau d’indice, si oui alors on le retourne, sinon on va dérouler la transformée inverse jusqu’à tombe sur un indice présent dans le tableau d’indice et on retourne la position correspondante plus le nombre d’étapes de transformées inverses effectués. L’algorithme est le suivant :\n",
+ "\n",
+ "\n",
+ "\n",
+ "ou k est l’indice recherché, idxs le tableau d’indices et stepleft est la fonction suivante :\n",
+ "\n",
+ "\n",
+ "\n",
+ "Notez que le modulo apparait car on peut se retrouver à dépasser la taille de la transformée initiale dans la recherche de l’indice le plus proche.\n",
+ "\n",
+ "Ecrivez la fonction index_from_match(k, idxs) -> int qui retourne l'indice dans la séquence d'origine à partir de k l'indice dans la transformée et idxs le tableau des correspondances des positions."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {
+ "id": "3RzUM2FJ91u_"
+ },
+ "outputs": [],
+ "source": [
+ "print(\"Votre code ici !\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "PtBlbmeNtaA2"
+ },
+ "source": [
+ "## Exercice 5 : Transformée progressive\n",
+ "\n",
+ "Le calcul de la transformée BW a un énorme coût en mémoire (il faut déterminer l’ordre des préfixes pour chaque position) ce qui la rend inutilisable pour des génomes. Une version légèrement modifiée permet d’échanger de la consommation mémoire contre du temps de calcul en séparant l’ensemble des suffixes à trier en k sous-ensembles, générer la transformée pour chaque sous-ensemble et les concaténer. Pour générer les sous-ensembles, on va simplement tirer k positions aléatoires le long de la séquence.\n",
+ "\n",
+ "Ecrivez la fonction BWT_prog(s, k) -> str qui retourne la transformée de la séquence s en utilisant k graînes.\n"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {
+ "id": "aoWF2thh8zEn"
+ },
+ "outputs": [],
+ "source": [
+ "print(\"Votre code ici !\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {
+ "id": "RsY4KN-x83cM"
+ },
+ "source": [
+ "## Exercice 6: reconstruction de génome\n",
+ "\n",
+ "Téléchargez le génome de la mitochondrie humaine: https://www.dropbox.com/scl/fi/ux5m9r80ejm4lz7w3qv0c/sequence_mito.fasta?rlkey=qoaz38vwmg3fgkymx09c7m59j&dl=0\n",
+ "Téléchargez l'ensemble des reads ici: https://www.dropbox.com/scl/fi/uajxylv3yfo1is876x1e9/reads_sequence_mito_l69_seed43_N2000_mu0_norc.fasta?rlkey=h16uoh8xjostqbbmy1tc1k86d&dl=0\n",
+ "\n",
+ "1. Calculez la transformée BW du génome\n",
+ "2. Alignez chaque read sur le génome\n",
+ "3. Calculez la séquence concessus"
+ ]
+ }
+ ],
+ "metadata": {
+ "colab": {
+ "authorship_tag": "ABX9TyO0SUJI6kaczFcOh8NoKcqb",
+ "include_colab_link": true,
+ "provenance": []
+ },
+ "kernelspec": {
+ "display_name": "Python 3 (ipykernel)",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.10.12"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}