Probabilité de trouver une séquence particulière de paires de bases

10

Penser à la probabilité me fait toujours réaliser à quel point je suis mal à compter ...

Considérons une séquence de lettres de base A ,n , chacun également susceptibles d'apparaître. Quelle est la probabilité que cette séquence contienne une séquence particulière de paires de bases d'intérêt de longueur r n ?A,T,C, and Grn

Il existe séquences différentes (tout aussi probables) possibles. Commencez par la séquence d'intérêt au début de la séquence complète; 4 n - r séquences comme celle-ci sont possibles. Nous pouvons commencer notre séquence d'intérêt dans n + 1 - r emplacements différents. Par conséquent, ma réponse est ( n + 1 - r ) / 4 r .4n4nrn+1r(n+1r)/4r

Cette probabilité augmente en , ce qui me semble logique. Mais cette probabilité dépasse 1 lorsque n > 4 r + r - 1 . Mais ça ne peut pas être. La probabilité devrait approcher 1 dans la limite (me semble), mais pas la dépasser.nn>4r+r1

Je suppose que je compte deux fois quelque chose. Qu'est-ce que je rate? Merci.

(Pour info, pas des devoirs, juste un exemple de jouet en préparation aux examens. Une question posée par mon ami biologiste moléculaire.)

Charlie
la source
C'est exact, cela ne devrait pas dépasser un car cela violerait les axiomes de probabilité: books.google.com/…
Chris Simokat
1
(Vaguement) lié: stats.stackexchange.com/questions/12174/…
cardinal

Réponses:

5

Considérons une petite version de ce problème avec . Quelle est la probabilité qu'une séquence de cinq lettres contienne la cible A C G T ? C'est facile: 4 - 4 de toutes les séquences commencent par cette chaîne, un autre 4 - 4 fin avec elle, et aucune séquence à la fois commence et se termine par cette chaîne. Par conséquent, la chance est de 2 × 4 - 4 .n=5ACGT44442×44

D'un autre côté, quelle est la chance de ? Encore une fois, 4 - 4 des séquences commencent par cette chaîne, la même fin de proportion avec cette chaîne, et 4 - 5 de toutes les séquences font les deux . Par conséquent, selon le principe d'inclusion-exclusion, la réponse est 2 × 4 - 4 - 4 - 5 .AAAA44452×4445

En général, la réponse dépend de la structure de la sous-chaîne. Pour être plus précis, lorsque vous numérisez une chaîne (de gauche à droite, par exemple) pour , vous ignorez tous les caractères jusqu'à ce que vous voyiez ce A initial . Après cela, il y a trois possibilités: le caractère suivant est une correspondance pour C , le suivant est une non-correspondance pour C mais n'est pas un A (vous êtes donc de retour dans l'état d'attente pour un A ), ou le suivant est un non-match mais c'est un A , vous plaçant dans l'état juste-vu- un - A . En revanche, envisagez une recherche pour A C T A C GACGTACCAAAAACTACG. Supposons que vous avez vu le préfixe . Le caractère suivant correspondra si elle est G . Quand c'est un non-match, (i) un C vous met dans l'état d'attente initial pour A , (ii) un A vous fait surveiller un C , et (iii) un T signifie que vous avez déjà vu A C T et vous êtes déjà à mi-chemin d'un match (et vous cherchez le deuxième A ). La "structure" pertinente consiste évidemment en des motifs de sous-chaînes dans la cible qui correspondent au préfixe de la cible. C'est pourquoi les chances dépendent de la chaîne cible.ACTACGCAACTACTA

Les diagrammes de FSA que je préconise dans une réponse à Time pris pour frapper un motif de têtes et de queues dans une série de lancers de pièces peuvent aider à comprendre ce phénomène.

whuber
la source
3

Une approximation grossière serait . Vous prenez la probabilité que votre séquence ne se produise pas à un endroit particulier, la mettez à la puissance du nombre d'emplacements (supposant faussement l'indépendance), qui est n - r + 1 et non n - r , et ceci est une approximation de son ne se produit pas, vous devez ensuite soustraire cela de 1 . 1(11/4r)nr+1nr+1nr1

Un calcul précis dépendra du motif précis que vous recherchez. est plus susceptible de se produire pas que A T C G T .AAAAAATCGT

Henri
la source
Peut - être juste moi, mais semble un peu plus clair en termes de comprendre comment l'équation a été construite. 1(1(1/4)r)n(r1)
@JoeRocc - Je soupçonne que c'est personnel. Si vous lisez de la page à la page 400 d'un livre, avez-vous lu 400 - 300 + 1 = 101 pages ou 400 - ( 300 - 1 ) = 101 pages? 300400400300+1=101400(3001)=101
Henry
Pas de soucis, je n'allais que par mon intuition du problème. Si nous dérivons intuitivement une équation pour être , alors en essayant de l'expliquer à quelqu'un, je pense qu'il vaut mieux le laisser comme ça plutôt que de le simplifier en a - b + c - 1 + d (bien que cela puisse certainement s'avérer plus intuitif après examen). Votre intuition peut avoir été différente dans tous les cas :)(une-(b-(c-1+)))une-b+c-1+
2

Vous comptez deux fois les séquences qui incluent plusieurs fois votre sous-séquence cible, par exemple à la fois à la position A et à la position B! = A. C'est pourquoi votre probabilité erronée peut dépasser 1

user145136
la source
Très bien fait ! +1
Michael R. Chernick
1

Il est possible d'obtenir la probabilité exacte d'une sous-séquence particulière en utilisant une représentation en chaîne de Markov du problème. Les détails de la façon de construire la chaîne dépendent de la sous-séquence particulière d'intérêt, mais je vais donner quelques exemples de la façon de procéder.


Probabilité exacte via la chaîne de Markov: considérons une séquence discrète de résultats de UNE,T,C,g où les résultats de la séquence sont échangeables, et supposons que nous soyons intéressés par une sous-chaîne de longueur k . Pour toute valeur donnée de n , disons W être le cas où la sous - chaîne d'intérêt se produit, et soit Hune être le cas où ces dernières une des résultats sont les premiers une<k caractères de la chaîne d'intérêt (mais pas plus que cela) . Nous utilisons ces événements pour donner la partition suivante de k+1 états d'intérêt possibles:

État 0W¯H0,   État 1W¯H1,   État 2W¯H2,   État 3W¯H3,   Etat k-1W¯Hk-1,Etat kW.  

Puisque la séquence des résultats est supposée être échangeable, nous avons des résultats indépendants conditionnels à leurs probabilités respectives θUNE+θT+θC+θg=1 . Votre processus d'intérêt peut être représenté comme une chaîne de Markov à temps discret qui commence dans l' État 0 à n=0 et transite selon une matrice de probabilité qui dépend de la sous-chaîne particulière d'intérêt. La matrice de transition sera toujours a (k+1)×(k+1)matrice représentant les probabilités de transition en utilisant les états ci-dessus. Si la sous-chaîne d'intérêt n'a pas été atteinte, chaque transition peut vous rapprocher de la sous-chaîne ou vous ramener à un état précédent qui dépend de la sous-chaîne particulière. Une fois la sous-chaîne atteinte, il s'agit d'un état absorbant de la chaîne, représentant le fait que l'événement d'intérêt s'est produit.

Par exemple, si la sous-chaîne d'intérêt est UNEUNEUNEUNEUNEUNE alors la matrice de transition est:

P=[1-θUNEθUNE000001-θUNE0θUNE00001-θUNE00θUNE0001-θUNE000θUNE001-θUNE0000θUNE01-θUNE00000θUNE0000001.]

Au contraire, si la sous-chaîne d'intérêt est UNECTUNEgC alors la matrice de transition est:

P=[1-θUNEθUNE00001-θUNE-θCθUNEθC00001-θUNE-θTθUNE0θT0001-θUNE000θUNE001-θUNE-θC-θgθUNEθC00θg01-θUNE-θCθUNE0000θC0000001.]

nP(W|n)={Pn}0,kn<k


Rn

#Create function to give n-step transition matrix for n = 1...N
#We will use the example of the substring of interest "AAAAAA"

#a is the probability of A
#t is the probability of T
#c is the probability of C
#g is the probability of G
#N is the last value of n
PROB <- function(N,a,t,c,g) { TOT <- a+t+c+g;
                              a <- a/TOT; 
                              t <- t/TOT; 
                              c <- c/TOT; 
                              g <- g/TOT; 

                              P <- matrix(c(1-a, a, 0, 0, 0, 0, 0,
                                            1-a, 0, a, 0, 0, 0, 0,
                                            1-a, 0, 0, a, 0, 0, 0,
                                            1-a, 0, 0, 0, a, 0, 0,
                                            1-a, 0, 0, 0, 0, a, 0,
                                            1-a, 0, 0, 0, 0, 0, a,
                                              0, 0, 0, 0, 0, 0, 1),
                                          nrow = 7, ncol = 7, 
                                          byrow = TRUE);
                              PPP <- array(0, dim = c(7,7,N));
                              PPP[,,1] <- P;
                              for (n in 2:N) { PPP[,,n] <- PPP[,,n-1] %*% P; } 
                              PPP }

#Calculate probability for N = 100 for equiprobable outcomes
N <- 100;
a <- 1/4;
t <- 1/4;
c <- 1/4;
g <- 1/4;
PROB(N,a,t,c,g)[1,7,N];

[1] 0.01732435

UNEUNEUNEUNEUNEUNEn=1000,01732435

Ben - Réintègre Monica
la source