Discussion :

[Jasmine80]

Bonjour,


Voici mon problème :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
my $sql0 = <<"SQL0";
SELECT accession, organism, sequence
FROM ma_table
WHERE sequence REGEXP 'TGTAGACGTGACGATGTAACAC[ATCG]{100,}' OR sequence REGEXP '[ATCG]{100,}TGATACCGATCCCTAGAACAGAT' 
SQL2

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

my ($subseq) = $sequence =~ m/((?:TGTAGACGTGACGATGTAACAC)?[ATCG]{100,250}(?:TGATACCGATCCCTAGAACAGAT)?)/;

$sequence contient TGTAGACGTGACGATGTAACAC ou TGATACCGATCCCTAGAACAGAT ou les 2 polynucléotides séparés par 100 à 250 nucléotides.

J'aimerais récupérer des sous-séquences de $sequence dont les 2 polynucléotides (TGTAGACGTGACGATGTAACAC et TGATACCGATCCCTAGAACAGAT) seraient les bornes. Si un seul polynucléotide est présent on récupère ce qu'il est possible de récupérer en gardant les même bornes.

Le problème est que si je fais une recherche gourmande, je récupère d'office 250 nucléotides après le premier polynucléotide (ce qui peut dépasser le second polynucléotide) et dans le cas contraire ([ATCG]{100,250}?), je n'en récupère que 100 et le second polynucléotide est laissé de côté.

Comment puis-je faire en une seule regexp?

Merci pour votre aide,

[Philou67430]

Erreur

[Philou67430]

As-tu essayé quelque chose du genre :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

my ($subseq) = $sequence =~ m/((?:^|TGTAGACGTGACGATGTAACAC)[ATCG]{100,250}(?:TGATACCGATCCCTAGAACAGAT|$))/;

[s.n.a.f.u]

As-tu essayé quelque chose du genre :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

my ($subseq) = $sequence =~ m/((?:^|TGTAGACGTGACGATGTAACAC)[ATCG]{100,250}(?:TGATACCGATCCCTAGAACAGAT|$))/;

Je ne peux pas tester, mais je pense qu'il faut utiliser le caractère ungreedy pour la recherche des nucléotides.
En l'état, sur une chaîne de longueur maximum, le motif va trouver 250 nucléotides, qu'il y ait présence ou non du polynucléotide à l'intérieur. Ensuite il va éventuellement backtracker ( beurk, c'est laid comme mot) et trouvera le polynucléotide de fin, mais pas celui qui est à l'intérieur.
Bref, je m'explique comme un manche, mais la version ungreedy est à tester :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

my ($subseq) = $sequence =~ m/((?:^|TGTAGACGTGACGATGTAACAC)[ATCG]{100,250}?(?:TGATACCGATCCCTAGAACAGAT|$))/;

Jasmine, tu avais bien cerné le problème de gourmandise (normal, tu es une femme... ), et ton souci venait principalement du caractère optionnel du dernier polynucléotide. Enfin je pense, j'ai pas matière à tester...

[Jasmine80]

J'avais résolu le problème en utilisant

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
	$sequence =~ s/.*(TGTAGACGTGACGATGTAACAC[ATCG]{100,250}).*/$1/;
	$sequence =~ s/.*?([ATCG]{100,250}TGATACCGATCCCTAGAACAGAT).*/$1/;

... est-ce moins performant?

Je vais tester vos réponses, merci pour votre aide.

[s.n.a.f.u]

Maintenant que j'y regarde une deuxième fois, je vois un autre problème potentiel :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(?:^|TGTAGACGTGACGATGTAACAC)

En début de chaîne, il va forcément trouver ^ et ne cherchera pas le polynucléotide.

[s.n.a.f.u]

Tu disais

$sequence contient TGTAGACGTGACGATGTAACAC ou TGATACCGATCCCTAGAACAGAT ou les 2 polynucléotides séparés par 100 à 250 nucléotides.

mais vue ta solution, il semble que cette séquence n'est qu'une partie de la chaîne, c'est bien cela ?

[Philou67430]

Tu as bien raison snafu, quand on n'a pas matière à tester, on arrive à écrire n'importe quoi

[Jasmine80]

Code s.n.a.f.u;4471190 :

Sélectionner tout - Visualiser dans une fenêtre à part

(?:^|TGTAGACGTGACGATGTAACAC)

Oui, il suffit d'inverser

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(?:GCATATCAATAAGCGGAGGA|^)

Jasmine, tu avais bien cerné le problème de gourmandise (normal, tu es une femme... ),

surtout le chocolat

mais vue ta solution, il semble que cette séquence n'est qu'une partie de la chaîne, c'est bien cela ?

Oui, normalement, c'est bien cela mais tous les cas sont possibles avec des chaînes tronquées.

[Jasmine80]

Le problème est que parfois on retrouve TGTAGACGTGACGATGTAACAC suivit de plus de 250 nucléotides sans TGATACCGATCCCTAGAACAGAT et donc cela pose problème. Il faudrait qu'il s'arrête au 250ième si iol n'y a pas TGATACCGATCCCTAGAACAGAT et si on n'est pas à la fin :

Peut-être ça :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

my ($subseq) = $sequence =~ m/((?:TGTAGACGTGACGATGTAACAC|^)[ATCG]{100,250}?(?:TGATACCGATCCCTAGAACAGAT|\w))/;

.. ça récupère le caractère (\w) en trop mais c'est un moindre mal.

[s.n.a.f.u]

Et que donnerai ceci :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
 
my ($subseq) = $sequence =~ m/^
.*?                               # élimination des inopportuns en début de chaîne
(                                 # début de capture
(?: TGTAGACGTGACGATGTAACAC)?      # recherche du premier polynucleotide facultatif
( (?: [ATCG](?!\2) ){100,250})   # parcours des nucléotides intermédiaires
(TGATACCGATCCCTAGAACAGAT)?        # deuxième polynucléotide facultatif
)                                 # fin de capture
/x;

[Jasmine80]

Et que donnerai ceci :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
 
my ($subseq) = $sequence =~ m/^
.*?                               # élimination des inopportuns en début de chaîne
(                                 # début de capture
(?: TGTAGACGTGACGATGTAACAC)?      # recherche du premier polynucleotide facultatif
( (?: [A-Z](?!=\2) ){2,100})      # parcours des nucléotides intermédiaires
(TGATACCGATCCCTAGAACAGAT)?        # deuxième polynucléotide facultatif
)                                 # fin de capture
/x;

Que signifie la ligne?

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

( (?: [A-Z](?!=\2) ){2,100})

Les séquences sont trop courtes ... J'ai un RDV à 16h30, je dois quitter mais je m'y remets dès demain matin.

Merci beaucoup pour ton aide.

[s.n.a.f.u]

Que signifie la ligne?

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

( (?: [A-Z](?!=\2) ){2,100})

Les séquences sont trop courtes ... J'ai un RDV à 16h30, je dois quitter mais je m'y remets dès demain matin.

Merci beaucoup pour ton aide.

J'ai corrigé, c'était mon cas de test parce que je n'ai pas beaucoup de nucléotides dans le bureau...

[Philou67430]

Une piste : peut-être est-il possible d'utiliser les expressions conditionnelles pour rendre le nucléotide de fin obligatoire si le nucléotide de début est vide...

(?(1)(?:TGATACCGATCCCTAGAACAGAT)?|TGATACCGATCCCTAGAACAGAT)

[s.n.a.f.u]

Je m'emmêle dans les conditionnels et autres facultatifs !

Finalement, ne serait-ce pas plus clair de gérer ça au cas par cas :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
 
#!/usr/bin/perl
use warnings;
use Data::Dumper;
 
my $start = 'AA';
my $end   = 'GG';
my $min = 2;
my $max = 5;
 
my $reg = "
(
    (?: ($start) ([ACGT]{$min,$max})? ($end) )
    |
    (?: ($start) ([ACGT]{$min,$max}) )
    |
    (?: ([ACGT]{$min,$max}) ($end) )
    |
    (?: ($start | $end) )
)";
 
 
$str = 'ACGAAGGGAC';
match($str, ++$try);
 
$str = 'ACGAAACGTAGGGAC';
match($str, ++$try);
 
$str = 'ACGAAACGTAGAC';
match($str, ++$try);
 
$str = 'ACGACGTAGGGAC';
match($str, ++$try);
 
$str = 'ACGAAGAC';
match($str, ++$try);
 
$str = 'ACGGGGAC';
match($str, ++$try);
 
$str = 'ACGAAACACACGTAGGGAC';
match($str, ++$try);
 
$str = 'ACGACACACGTAGGGAC';
match($str, ++$try);
 
sub match {
    my ($str, $try) = @_;
    print "Essai n° $try\n";
    @res = $str =~ m/$reg/x;
    print Dumper(@res);
    print "\n";
}

J'ai changé les valeurs des (poly)nucléotides pour tests, mais il suffit de paramétrer correctement pour tes chaînes.
A suivre...

[Jasmine80]

Je m'emmêle dans les conditionnels et autres facultatifs !

Finalement, ne serait-ce pas plus clair de gérer ça au cas par cas

Oui, ça devient un vrai casse-tête. Je vais garder mes 2 substitutions :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
	$sequence =~ s/.*(TGTAGACGTGACGATGTAACAC[ATCG]{100,250}).*/$1/;
	$sequence =~ s/.*?([ATCG]{100,250}TGATACCGATCCCTAGAACAGAT).*/$1/;

... le programme prend quelques secondes sur mes 638 séquences ... je pensais qu'il y avait plus simple et plus performant en une seule regexp mais apparemment pas Je préfère un programme un peu moins performant mais compréhensible.


Que signifie

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(?: [ATCG](?!\2) ){100,250}

... je ne comprends pas le (?!\2).

Une piste : peut-être est-il possible d'utiliser les expressions conditionnelles pour rendre le nucléotide de fin obligatoire si le nucléotide de début est vide...

(?(1)(?:TGATACCGATCCCTAGAACAGAT)?|TGATACCGATCCCTAGAACAGAT)

A quoi sert le (1)?


Merci pour votre aide à tous les 2.

[Philou67430]

je ne comprends pas le (?!\2).

c'est l'assertion négative de taille zero de recherche par avance (traduction approximative).

"(?!pattern)"
A zero-width negative look-ahead assertion. For example
"/foo(?!bar)/" matches any occurrence of "foo" that isn't
followed by "bar". Note however that look-ahead and
look-behind are NOT the same thing. You cannot use this
for look-behind.

If you are looking for a "bar" that isn't preceded by a
"foo", "/(?!foo)bar/" will not do what you want. That's
because the "(?!foo)" is just saying that the next thing
cannot be "foo"--and it's not, it's a "bar", so "foobar"
will match. You would have to do something like
"/(?!foo)...bar/" for that. We say "like" because
there's the case of your "bar" not having three
characters before it. You could cover that this way:
"/(?:(?!foo)...|^.{0,2})bar/". Sometimes it's still
easier just to say:

if (/bar/ && $` !~ /foo$/)

For look-behind see below.

A quoi sert le (1)?

C'est plus simple, c'est le pattern conditionnel, en l'occurence, conditionné par le match de la première capture du motif.

"(?(condition)yes-pattern|no-pattern)"
"(?(condition)yes-pattern)"
Conditional expression. "(condition)" should be either an
integer in parentheses (which is valid if the corresponding
pair of parentheses matched), a
look-ahead/look-behind/evaluate zero-width assertion, a name
in angle brackets or single quotes (which is valid if a
buffer with the given name matched), or the special symbol
(R) (true when evaluated inside of recursion or eval).
Additionally the R may be followed by a number, (which will
be true when evaluated when recursing inside of the
appropriate group), or by &NAME, in which case it will be
true only when evaluated during recursion in the named group

Here's a summary of the possible predicates:

(1) (2) ...
Checks if the numbered capturing buffer has matched
something.
...

Merci pour votre aide à tous les 2.

Si une expression régulière unique existe, il n'est pas évident qu'il puisse être plus efficace à cause du backtracking...

[s.n.a.f.u]

Que signifie

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(?: [ATCG](?!\2) ){100,250}

... je ne comprends pas le (?!\2).

Salut,

en fait, il y a deux concepts dans ce hiéroglyphe :

(?! expression ) : c'est une "assertion négative avant", ce qui peut se traduire plus simplement par "qui n'est pas suivi de expression".

\2 : c'est tout simplement le contenu capturé par la deuxième paire de parenthèses

Voilou !

[s.n.a.f.u]

Oui, ça devient un vrai casse-tête. Je vais garder mes 2 substitutions :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
	$sequence =~ s/.*(TGTAGACGTGACGATGTAACAC[ATCG]{100,250}).*/$1/;
	$sequence =~ s/.*?([ATCG]{100,250}TGATACCGATCCCTAGAACAGAT).*/$1/;

Alors là je dirais oui, mais non.

Car avec cette substitution, tu fais une hypothèse de plus : les 100 à 250 nucléotides intermédiaires sont obligés d'être présents, ainsi que les deux bornes.
Du coup, ça simplifie énormément le problème, mais ça ne prend pas en compte tous les cas de départ.
Ca ne répond par exemple pas au cas où les deux bornes sont collées l'une à l'autre sans rien entre les deux.

[Jasmine80]

Merci beaucoup Philou pour ces explications j'apprends plein de choses sur ce forum

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
my ($subseq) = $sequence =~ m/^
.*?                               # élimination des inopportuns en début de chaîne
(                                 # début de capture
(?: TGTAGACGTGACGATGTAACAC)?      # recherche du premier polynucleotide facultatif
( (?: [A-Z](?!=\2) ){100,250})      # parcours des nucléotides intermédiaires
(TGATACCGATCCCTAGAACAGAT)?        # deuxième polynucléotide facultatif
)                                 # fin de capture
/x;

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(?: [A-Z](?!=\2) ){100,250}

On recherche donc une majuscule qui n'est pas suivie par =\2 est-ce bien le contenu de $2 ... j'ai du mal à comprendre le =.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(?(1)(?:TGATACCGATCCCTAGAACAGAT)?|TGATACCGATCCCTAGAACAGAT)

Donc si on matches le premier TGATACCGATCCCTAGAACAGAT, (1) est valide et on ne recherche pas le second motif ... qui devrait-être le second polynucléotide non?