Aujourd'hui, nous allons recoder la fonction strstr.
Mise en situation
La recherche en génétique progresse chaque jour. Elle a déjà donné d'excellents résultats, et continue d'être active. Seulement, la quantité d'informations contenue dans un brin d'ADN est simplement énorme... Il est donc logique de la traiter par ordinateur !
Les séquences ADN peuvent être représentées par un alphabet de quatre lettres : A, T, G ,C. Pour connaître la signification de ces lettres, je vous invite à consulter la page wikipédia concernant l'ADN. Une opération courante dans le domaine est la recherche d'un gène (qui sera ici une séquence de lettres) dans une séquence ADN bien plus grande.
Nous allons donc coder un programme capable de retrouver une sous-chaîne dans un texte !
Énoncé
Écrivez un programme qui lit un gène sur la première ligne, une séquence ADN sur la deuxième, puis qui donne l'indice ou démarre le gène dans la séquence.
Exemple :
ATGC
ATTGGCCTATGCGGA
En sortie :
9
Niveau 1
Écrivez le programme de recherche en utilisant l'algorithme le plus simple. Bien entendu, il est interdit de recourir à la fonction strstr, le but est que vous la fassiez vous-mêmes.
Niveau 2
Améliorez votre programme de recherche pour qu'il utilise l'algorithme de Rabin-Karp. L'idée est la suivante : nous allons représenter chaque lettre par un chiffre. Par exemple, A=0, T=1, G=2 et C=3. Nous pouvons alors écrire un gène, tel que ATGC, comme un nombre entier, qui vaut 123.
Pour chercher le gène, il suffit de prendre les 4 premiers caractères de la séquence, d'opérer la même conversion, puis de comparer les résultats. Il nous faut ensuite convertir les caractères 2 à 5 de la séquence, et comparer à nouveau... Mais, lorsque nous avons converti les caractères 1 à 4 tout à l'heure, nous avons en particulier converti les caractères 2 à 4 en entiers ! N'y a-t-il pas moyen d'en profiter ?
Reprenons l'exemple précédent pour nous fixer les idées. Au premier tour, de boucle, la situation était :
0123
0112GCCTATGCGGA
Comme 123 et 112 sont différents, le gène n'est pas en position 1. Nous nous décalons donc d'une case dans la séquence :
_0123
A1122CCTATGCGGA
Donc, sachant que la séquence actuelle vaut 0112, et que le prochain caractère sera associé au nombre « 2 », comment en déduire que le prochain entier à comparer sera 1122 ?
Niveau 3
Le code précédent constitue une bonne amélioration de la recherche naïve codée au niveau 1. Cependant, nous sommes limités par la taille des entiers que nous pouvons manipuler ! En effet, comme vous le savez sûrement, une variable de type int ou long comporte une limite de grandeur pour les entiers stockés. Cette limite maximale varie selon la machine ; vous pouvez l'afficher grâce au code suivant.
#include <limits.h>
int main(void)
{
printf("%ul\n", ULONG_MAX);
return 0;
}
Ainsi, nous ne pouvons pas rechercher des gènes arbitrairement longs ! Il faut que leur « version numérique » tienne dans un long (ou autre type que vous aurez choisi). Nous avons deux solutions pour pallier à ce problème.
La première est d'utiliser un type de données capables de manipuler des nombres arbitrairement grands. L'objet de l'un des exercices du forum est d'écrire un tel type, vous pouvez donc le réutiliser si ça vous amuse.
Cependant, les opérations sur de grands nombres sont relativement lentes, et nous perdrions le bonus en rapidité que nous fait gagner l'algorithme de Rabin-Karp.
La seconde solution est d'utiliser un hash. L'idée est la suivante : nous allons appliquer un modulo au gène et au bloc à tester. Notez que ceci ne modifie en rien l'astuce du niveau 2 : elle fonctionne toujours.
En revanche, lorsque nous trouverons deux nombres égaux, ne ne seront plus complètement sûrs que les gènes correspondants le sont bien aussi (du fait du modulo). Il devient donc nécessaire de faire une comparaison textuelle en cas d'égalité. Par contre, si les deux nombres trouvés sont différents, les deux séquences sont nécessairement différentes, et il n'y a plus besoin de test (ce qui est un peu l'intérêt quand même).
Conclusion
À l'issue de cet exercice, nous avons implémenté un algorithme efficace de recherche de sous-chaîne, qui évite « souvent » les comparaisons inutiles.
Il existe encore de nombreuses autres méthodes pour pour rechercher efficacement une sous-chaîne dans une chaîne. Vous pensiez que l'exercice de M@teo21 était fini, enterré, oublié ? Vous avez pourtant de quoi vous amuser encore longtemps...
Si vous souhaitez plus de renseignements sur les défis, rendez vous sur le topic de recensement des défis, vous y trouverez également les règles principales.
Cet exercice a été écrit par GuilOooo, merci à lui !
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
int recherche(char* motif, char* chaine){
int i, k=0;
int tailleMotif = strlen(motif);
int tailleChaine = strlen(chaine);
for(i=0; i<(tailleChaine-tailleMotif+1); i++){
for(k=0; (k<tailleMotif)&&(chaine[i+k] == motif[k]); k++);
if (k==tailleMotif)
return (i+k)-tailleMotif+1;
}
return -1;
}
int main(int argc, char** argv)
{
char* motif = malloc(sizeof(char)*50);
char* chaine = malloc(sizeof(char)*50000);
printf("Entrez le motif : ");
scanf("%s", motif);
printf("Entrez la chaine : ");
scanf("%s", chaine);
printf("Resultat : %d", recherche(motif, chaine));
free(motif);
free(chaine);
return 0;
}
Niveau 2&3
Le code est assez moche, mais je me demande si on peut améliorer le calcul à la ligne 91.
Et je suis pas sur d'avoir mis le modulo au bon endroit
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <limits.h>
#define TRUE 1
#define FALSE 0
#define MODULO INT_MAX
int recherche(char* motif, char* chaine){
int i, k=0;
int tailleMotif = strlen(motif);
int tailleChaine = strlen(chaine);
for(i=0; i<(tailleChaine-tailleMotif+1); i++){
for(k=0; (k<tailleMotif)&&(chaine[i+k] == motif[k]); k++);
if (k==tailleMotif)
return (i+k)-tailleMotif+1;
}
return -1;
}
int convert(char carac){
//Transforme un caractere de nucléotide dans son entier correspondant
int lettre;
switch (carac){
case 'A': lettre = 0; break;
case 'T': lettre = 1; break;
case 'G': lettre = 2; break;
case 'C': lettre = 3; break;
default: lettre = 0; break;
}
return lettre;
}
int hash(char* str){
//Transforme une chaine de nucléotide en son entier correspondant
int i;
int tailleStr = strlen(str);
int lettre, resultat = 0;
for(i=0; i<tailleStr; i++){
lettre=convert(str[i]);
resultat += lettre*puiss10(tailleStr-i-1);
}
return resultat%MODULO;
}
int strCmp(char* str1, char* str2){
//Histoire de pas trop utiliser string.h non plus ...
int i = 0;
if (strlen(str1) == strlen(str2)){
for(i=0; i<strlen(str1); i++){
if(str1[i] != str2[i])
return FALSE;
}
return TRUE;
}
return FALSE;
}
int puiss10(int puissance){
//Retourne 10^puissance
int resultat = 1;
for(; puissance>0;puissance--)
resultat *= 10;
return resultat;
}
int rechercheRabinKarp(char* motif, char* chaine){
int i;
int tailleMotif = strlen(motif);
int tailleChaine = strlen(chaine);
int hashMotif = hash(motif);
//On extrait une sous-chaine du début de la séquence,
//de la taille du motif
char* sousChaine = malloc(tailleMotif+1);
strncpy(sousChaine, chaine, tailleMotif);
sousChaine[tailleMotif] = '\0';
//Fin de l'extraction
int hashSousChaine = hash(sousChaine);
for(i=0; i<(tailleChaine-tailleMotif+1); i++){
if(hashMotif == hashSousChaine){
strncpy(sousChaine, chaine+i, tailleMotif);
sousChaine[tailleMotif] = '\0';
if(strCmp(motif, sousChaine))
return i+1;
}
hashSousChaine = ((hashSousChaine%(puiss10(tailleMotif-1))*10
+
convert(chaine[i+tailleMotif])))%MODULO;
}
free(sousChaine);
return -1;
}
int main(int argc, char** argv)
{
char* motif = malloc(50000);
char* chaine = malloc(500000);
printf("Entrez le motif : ");
scanf("%s", motif);
printf("Entrez la chaine : ");
scanf("%s", chaine);
printf("Resultat : %d", rechercheRabinKarp(motif, chaine));
free(motif);
free(chaine);
return 0;
}
C'est moi ou il n'y a que peu de différence entre les niveau 2 et 3(à un modulo et un strncmp près)? Le modulo de plus s'effectue naturellement sur des unsigned.
Pour le niveau 3
#include <stdio.h>
#include <string.h>
struct t_hash
{
unsigned long value;
unsigned long mod;
};
static unsigned long powOf10(unsigned long n)
{
unsigned long ret = 1;
for(unsigned long i = 0; i < n; i++)
ret *= 10;
return ret;
}
static size_t getIndex(char c)
{
switch(c)
{
case 'A':
return 1;
case 'C':
return 2;
case 'T':
return 3;
case 'G':
return 4;
default:
return 0;
}
}
struct t_hash hash_init(char const *s, size_t n)
{
struct t_hash ret = {0, 0};
ret.mod = powOf10(n - 1);
for(size_t i = 0; s[i] && i < n; i++)
ret.value = (ret.value * 10) + getIndex(s[i]);
return ret;
}
void hash_update(struct t_hash *hash, char c)
{
hash->value = ((hash->value % hash->mod) * 10) + getIndex(c);
}
char *find(char const *src, char const *needle)
{
size_t len = strlen(needle);
struct t_hash hash_needle = hash_init(needle, len);
struct t_hash hash_crt = hash_init(src, len);
for(char const *nxt_char = src + len - 1; *nxt_char; src++)
{
if(hash_needle.value == hash_crt.value && strncmp(src, needle, len) == 0)
return (char *)src;
hash_update(&hash_crt, *++nxt_char);
}
return NULL;
}
int main(void)
{
char const *src[5] = {"", "ATGC", "AATGC", "AAATGC", "ATG"};
char const *needle = "ATGC";
for(size_t i = 0; i < 5; i++)
{
char *p = find(src[i], needle);
if(p)
printf("\"%s\" dans \"%s\" en position %zu\n", needle, src[i], p - src[i]);
else
printf("\"%s\" pas dans \"%s\"!\n", needle, src[i]);
}
return 0;
}
Ce serait pas mal de fournir un fichier de test, histoire qu'on fasse des tests communs.
Tiens intéressant, il y a du monde sur ce défis, allez un petit tour des codes présenté
Bon un petit tour des codes présenté
@lucas-84: Humm :
ATT
TAATTF
/
@redspikers: Pas de problème niveau implémentation, l'exercice niveau 1 est correcte Je regarderais le niveau 3 plus tard ^^"
Après, je me permet quand même de faire 2-3 commentaires rapide :
L'allocation dynamique n'est pas franchement nécessaire ici ... Surtout si tu ne vérifies pas les retours de malloc
sizeof(char) est toujours égal à 1 donc sizeof(char) * i est toujours égal à i
Tu pourrais peut-être éviter de chercher la taille de la chaine avec strlen, ce sont des parcours inutiles
@redspikers: Pas de problème niveau implémentation, l'exercice niveau 1 est correcte Je regarderais le niveau 3 plus tard ^^"
Après, je me permet quand même de faire 2-3 commentaires rapide :
L'allocation dynamique n'est pas franchement nécessaire ici ... Surtout si tu ne vérifies pas les retours de malloc
Pour l'allocation dynamique, ça me permet d'avoir des tableaux plus grand. Exemples :
Chez moi : </ul>
char chaine2[5000000]; //plante alors que :
char* chaine2 = malloc(5000000); //ne plante pas.
Et pourtant je demande la même taille.
sizeof(char) est toujours égal à 1 donc sizeof(char) * i est toujours égal à i</li>
Yep, je l'ai fais par réflexe, mais je l'ai pas fait pour le niveau 2&3
Tu pourrais peut-être éviter de chercher la taille de la chaine avec strlen, ce sont des parcours inutiles
</li> La je vois pas pourquoi c'est inutile ? D'autant plus que strlen je le fais que au début de la recherche.for(k=0;(k<tailleMotif)&&(chaine[i+k]==motif[k]);k++);// <3 </li>
Je pensais que j'allais me faire insulter pour avoir trop voulu compacter
char chaine2[5000000]; //plante alors que :
char* chaine2 = malloc(5000000); //ne plante pas.
Normal On peut accéder à plus de mémoire sur le tas que par la pile.
Mais après, on pas besoin de 50 000 cases ^^" 200 c'est largement suffisamment ...
Citation : redspikers
Je pensais que j'allais me faire insulter pour avoir trop voulu compacter
Je vous donne mon code (niveau 1) que je viens de finir, si vous avez des remarques à faire dessus n'hésites pas
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "lecture.h"
int main()
{
char gene[20] = {0};
char sequence[100] = {0};
int vrai = 0;
int i = 0, j = 0;
printf("Rentrez votre gene a rechercher : ");
lire(gene, 20);
printf("Rentrez la sequence de recherche : ");
lire(sequence, 1000);
while (vrai != strlen(gene))
{
vrai = 0;
for (i = 0; gene[i] != '\0'; i++)
{
if (gene[i] == sequence[j+i])
{
vrai+= 1;
}
}
j++;
}
printf("\nIl y a une occurence a la position %d", j);
return EXIT_SUCCESS;
}
(la fonction lire est celle étudiée dans le tutoriel de M@teo21, qui supprime le \n et vide le buffer).
En ce qui concerne le niveau 2 et 3 j'ai un peu de mal à comprendre le principe en fait : j'ai regardé un peu sur internet en quoi consistait une fonction de hachage mais j'ai pas vraiment compris Est-ce que dans ce cas-là on rentre "ATCG" par exemple et que la fonction de hachage va convertir ce ATCG en "0123" ou on rentre 0123 directement (je suppose que c'est la première solution) ?
Et sinon, là c'est juste à titre d'informations, à quoi cela sert de convertir cela en nombres pour comparer leurs valeurs numériques au lieu de directement comparer les caractères entre eux ? Si quelqu'un pourrait me répondre ce serait vraiment sympa
Et sinon, là c'est juste à titre d'informations, à quoi cela sert de convertir cela en nombres pour comparer leurs valeurs numériques au lieu de directement comparer les caractères entre eux ? Si quelqu'un pourrait me répondre ce serait vraiment sympa
Plus rapide à comparer. En plus, ça réduis le nombre de comparaisons. Si tu veux plus de détails, demande
Ton code ne donne pas une boucle infinie si le gène n'est pas trouvé ?
Aussi : lire(sequence,1000);
Avec : charsequence[100]={0};
Dernier point : while(vrai!=strlen(gene))
Niveau performance, c'est horrible ! En effet, on recalcule la taille de gène à chaque tour de boucle ...
🍊 - Étudiant - Codeur en C | Zeste de Savoir apprenez avec une communauté | Articles- ♡ Copying is an act of love.
Pour la boucle infinie oui, avec cette solution je n'en trouvais pas d'autre pour arrêter le programme s'il n'en trouvait pas parce que la séquence n'a pas réellement de taille maximale prédéfinie :/
Pour le lire(sequence, 1000) c'est juste un faute de frappe, j'ai mis un 0 en trop
Pour le while d'accord, je vais créer une variable qui stockera strlen(gene) dans ce cas.
Niveau 2 et 3, on peut optimiser un peu en utilisant la base 4 au lieu de la base 10.
Je ne sais pas si j'ai bien compris la technique du modulo, mais si oui, l’intérêt se retrouve limité dans le cas ou la séquence contient des gènes très ressemblants, et que l'on cherche un gène de longueur conséquente (donc demandant un modulo). En effet, le modulo permet seulement de savoir si la fin des 2 chaines (de longueur limitée, ici 31 caractères) sont équivalentes. Donc la fonction de hash n'en est pas véritablement une ici.
Mon code utilise les opérateurs bitwise à la place de l'addition et la multiplication. Le modulo est inutile, car le nombre est automatiquement tronqué par overflow.
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <limits.h>
long long unsigned val(char c)
{
switch(c)
{
case 'A':
return 0;
case 'T':
return 1;
case 'G':
return 2;
case 'C':
return 3;
}
fprintf(stderr, "Invalid sequence : %c found\n", c);
exit(1);
}
long long unsigned hash (const char* str, size_t len)
{
long long unsigned h = 0;
size_t i = 0;
for (i=0; i < len; ++i)
h = (h << 2) | val(str[i]);
return h;
}
long long unsigned secure_shift (long long unsigned n, size_t k)
{
return (k >= sizeof(long long unsigned) * CHAR_BIT ? 0 : n << k);
}
int karp_rabin(const char* seq, const char* gene)
{
size_t s_len = strlen(seq);
size_t g_len = strlen(gene);
if (s_len < g_len || g_len == 0)
return -1;
long long unsigned g_hash = hash(gene, g_len);
long long unsigned s_hash = hash(seq, g_len);
if (s_hash == g_hash && !strncmp(seq, gene, g_len))
return 0;
size_t i, j;
for (i=1, j=g_len; j < s_len; ++i, ++j)
{
s_hash = ( ( s_hash - secure_shift(val(seq[i-1]), 2*(g_len-1)) ) << 2 ) | val(seq[j]);
if (g_hash == s_hash && !strncmp(seq+i, gene, g_len))
return i;
}
return -1;
}
int main (void)
{
char* str = "ATTGGCCTATGCGGA";
char* pattern = "ATGC";
printf("%d\n", karp_rabin(str, pattern));
return 0;
}
EDIT : Il semble que l'exemple donné sur la page française de Wikipédia soit mauvais (induit en erreur), la page anglaise indique d'utiliser un nombre premier pour la base, ce qui semble plus logique.
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <limits.h>
#define BASE 5
long long unsigned val(char c)
{
switch(c)
{
case 'A':
return 0;
case 'T':
return 1;
case 'G':
return 2;
case 'C':
return 3;
}
fprintf(stderr, "Invalid sequence : %c found\n", c);
exit(1);
}
long long unsigned hash (const char* str, size_t len)
{
long long unsigned h = 0;
size_t i = 0;
for (i=0; i < len; ++i)
h = (h * BASE) + val(str[i]);
return h;
}
long long unsigned secure_pow (long long unsigned n, unsigned exp)
{
long long unsigned ret = 1;
while (--exp)
{
if (ret > ULLONG_MAX / BASE) // overflow
return 0;
ret *= n;
}
return ret;
}
int karp_rabin(const char* seq, const char* gene)
{
size_t s_len = strlen(seq);
size_t g_len = strlen(gene);
if (s_len < g_len || g_len == 0)
return -1;
long long unsigned g_hash = hash(gene, g_len);
long long unsigned s_hash = hash(seq, g_len);
if (s_hash == g_hash && !strncmp(seq, gene, g_len))
return 0;
long long unsigned p = secure_pow(BASE, g_len);
size_t i, j;
for (i=1, j=g_len; j < s_len; ++i, ++j)
{
s_hash = ((s_hash - (val(seq[i-1]) * p)) * BASE) + val(seq[j]);
if (g_hash == s_hash && !strncmp(seq+i, gene, g_len))
return i;
}
return -1;
}
int main (void)
{
char* str = "ATTGGCCTATGCGGA";
char* pattern = "ATGC";
printf("%d\n", karp_rabin(str, pattern));
return 0;
}
C'est gentil, GurneyH . Ton code à l'air sympa lui aussi (je suis assez minimaliste effectivement, mais ce n'est pas toujours une qualité...)
Au fait, j'ai fait des tests grandeur nature, mes deux codes sont faux. En fait, j'ai supposé à tort certains comportements d'overflow. J'ai corrigé le premier, reste à corriger le deuxième.
Ah, et au passage : dur dur de battre l'implémentation naïve dans les cas courants...
En fait, j'ai supposé à tort certains comportements d'overflow.
Tu peux préciser le problème stp? j'ai la flemme de lire les corrections apportée, j'avoue! Sur des unsigned on a bien un modulo en cas d'overflow, non?
Citation : yoch
(je suis assez minimaliste effectivement, mais ce n'est pas toujours une qualité...)
Bah minimaliste sûrement , mais j'apprends toujours beaucoup avec tes codes!
Tu peux préciser le problème stp? j'ai la flemme de lire les corrections apportée, j'avoue! Sur des unsigned on a bien un modulo en cas d'overflow, non?
En fait, j'ai travaillé sans utiliser de modulo, d’où mes soucis.
J'ai supposé que par exemple :
unsigned a = 1;
a << 100;
// ici, je pensais que a == 0, ce qui est faux...
Ensuite, pour mon second code (encore non corrigé), je pensais qu'un overflow de multiplication se comporte comme un modulo. Il me semble que ce n'est la cas qu'en base 2 (à confirmer que c'est bien le cas), mais pas dans les autres bases, donc on a bien besoin d'un modulo. (Au passage, même avec un modulo, il faut faire modulo une certaine puissance de la base pour avoir un résultat correct)
Plus sérieusement, je voulais vraiment passé par sscanf pour la conversion chaine -> nombre, d'où le code obtenu (et la nécessité de crée un format de manière dynamique)
× Après avoir cliqué sur "Répondre" vous serez invité à vous connecter pour que votre message soit publié.
× Attention, ce sujet est très ancien. Le déterrer n'est pas forcément approprié. Nous te conseillons de créer un nouveau sujet pour poser ta question.
🍊 - Étudiant - Codeur en C | Zeste de Savoir apprenez avec une communauté | Articles - ♡ Copying is an act of love.
🍊 - Étudiant - Codeur en C | Zeste de Savoir apprenez avec une communauté | Articles - ♡ Copying is an act of love.
🍊 - Étudiant - Codeur en C | Zeste de Savoir apprenez avec une communauté | Articles - ♡ Copying is an act of love.
🍊 - Étudiant - Codeur en C | Zeste de Savoir apprenez avec une communauté | Articles - ♡ Copying is an act of love.
🍊 - Étudiant - Codeur en C | Zeste de Savoir apprenez avec une communauté | Articles - ♡ Copying is an act of love.