DNA密码研究

Distribution of Codon Locations in DNA SequenceWu XiaoxiaCollege of Mathematics and Statistics, Wuhan University, Wuhan, 430072, Hubei, ChinaE-mail: wuxiaoxia_1015@ AbstractIn the analysis of DNA sequence, this paper introduces a method: stationary m-order Markov chain to study codon locations, from which the estimators of parameters can be obtained via its MLE. This paper gives a description of the exact distribution of the distance between codon occurrences along a DNA sequence and an asymptotic distribution of codon count through the Gaussian distribution. Keywords: codon; stationary m-order Markov chain; transition matrix. CLC number : O 2121 Introduction of the model in DNA sequenceA biological sequence is either a DNA sequence or a protein sequence, which is a finite sequence of letters either in the 4-letter DNA alphabet {},,,Αa c g t =. To model a biological sequence, we will consider models for random sequences of letters ()i i X X ∈= on a finite alphabet A , where is the set of integers.We consider a good model, the model Mm as an ergodic stationary m-order Markov chain on a finite alphabet A with transition matrix()()()()121121,,,,,m m m m a a a a a a a a Απ++∈Π=such that()()121111,|,,m m i m i m i m a a a a P X a X a X a π++−−====Define ()()12112,m i m i m a a a P X X a a a i µ+−==∀∈ , where µas a stationary distribution satisfy the equation µµ=Π, such that the equation()()()12121121,m m m m b a a a ba a a ba a a a Αµµπ−−∈=∑is satisfied for all 12()mma a a ∈ A .( see Ref.[1])A coding DNA sequence is naturally read as successive non-overlapping 3-letter words called codons. The codons are then translated into amino acids via the genetic code to produce a protein sequence .It is known that there are 64 codons for they are the compositions of three different letters from Α.The codons are used to be translated into amino acids except 3 terminators. And we know there are only20 amino acids. So there are several different codons coding for the same amino acid, and often the first two letters of some codon suffice to determine the corresponding amino acid. Therefore, letters may have different importance depending on their positions with respect to the codon partition. To distinguish the letter probabilities according to their positions modulo of 3 in the coding DNA sequence, we consider a stationary Markov chain with three distinct transition matrices 12,∏Πand3Π, such that for()121,m m a a a a Α+∈ and {}1,2,3k ∈,()()121313131,|,,,k m m j k m j k m j k m a a a a P X a X a X a j π++++−+−====∀∈This is model 3Mm . The index {}1,2,3k ∈is called phase and represents the position of aletter inside a codon. The stationary distribution k µon mΑ is given by()1231312(),k m j k m j k m a a a P X X a a a j µ+−++==∀∈such that the equations()(){}121121121(),,1,2,3k m k m k m m b a a a ba a a ba a a a k µµπ−−−∈Α=∈∑are satisfied for all 12()m m a a a Α∈ .Robin and Daudin [2] obtained the recursive formula for the exact distribution of the distance D between two word occurrences. It was first proposed for independent and uniformly distributed letters by Blom and Thorbum [3] . Similar results are derived in Robin and Daudin [4] and Stefanov [5] for the probability distribution of the distance between any word in a given set. The method of obtaining the exact distribution of the word count presented here generalizes the result that Gentleman and Mullin [6] obtained for case that the sequence is composed of i.i.d. letters, where each letter occurs with equal probability. The martingale and the conditional ones have advantage to provide explicit formulas for the asymptotic variance for the first order markov chain model. They are both due to Prum et al. [7] and to Schbath [8] for higher order models and Phased model. In this article, we consider the location of codon and apply the methods above to obtain its distribution.2 Estimation of the model parametersIn this paper we assume that 12nX X X is a stationary 1-order Markov chain on Αwithtransition matrix ()(),,a b a b π∈ΑΠ=and stationary distribution ()()a a µ∈Αfor simplicity.Generally, a DNA sequence contains a promoter before it starts to read the codons with hereditable information and ends with a terminator in the process of gene transcription. But in some other particular cases, it has not promoter, so the reading frame is not fixed and the results of transcription are distinct. Here we just consider the general case. That is it has a promoter and the reading frames are all fixed. So the sequence 12nX X X starts from the first codon and ends with aterminator. Here we assume that n can be divided by 3, otherwise we take []/33n n ′=×.Modeling a biological sequence, we need choose a probabilistic model and estimate the parameters according to the unique realization that is the biological sequence. We now give an estimator that maximizes the likelihood of 1M model. Assume 12nXX Xis a stationary Markovchain on A with transition matrix ()(),,a b a b π∈ΑΠ=and stationary distribution ()()a a µ∈Α. Thelikelihood L of the model is(){}()()()()()31,1,,,,1,2,3,k N ab k k a b k L a b a b k X a b πµπ∈Α=∈Α∈=∏∏where ()k N ab denotes the number of the 2-letter word ab which ends at301,1,2,33n i k i k +≤≤−=⎛⎞⎜⎟⎝⎠in the random sequence 12nX X X .To find the transition probabilities that maximize the likelihood, one maximizes the log likelihood(){}()()()()31,1log ,,,,1,2,3log log ,k kka b k L a b a b k X N ab a b ΑπΑµπ∈=∈∈=+∑∑We can separately maximize()()31log ,kkk b N ab a b π=∈Α∑∑for a Α∈, keeping in mind that(){},1,1,2,3kb a b k π∈Α=∀∈∑. Leta ∈Αand choose c Α∈, we have()()()()()()333111log ,log ,log 1,kkkkkkk b k b ck b cN ab a b N ab a b N ac a b Απππ=∈=≠=≠⎛⎞=+−⎜⎟⎝⎠∑∑∑∑∑∑ And for b c ≠,()()()()()()()31log ,,,,k k k k k b k k k N ab a b N ab N ac a b a b a c Αππππ=∈∂=−∂∑∑{},1,2,3k ∈.All the partial derivatives equal to zero means that ()()()(),,,k k k k N ab N ac b a b a c Αππ=∀∈. This implies inparticular that()()()()(),.,,kk d kd kkd N ad N ab N ad b a b a d Αππ∈Α∈Α∈Α==∀∈∑∑∑ It follows that ()()()()(){}1ˆ,,1,2,3.k k k kk d N ab N ab a b k N ad N a π−∈Α==∀∈⋅∑For notational convenience, theestimators mainly used in the reminder of the article will be ()()()1ˆ,k k k N ab a b N a π−= since ()()N a N a ⋅=except for the letter of the sequence for which the counts differ by 1.3 Exact distribution of the distance between codon occurrencesHere we are concerned with the length of the gaps between codon occurrences. Now wedescribe how to obtain the exact distribution of the distance between successive occurrences of a codon.Let 123ωωωω=be a word of length 3 on a finite alphabet A , where ωis a codon. We assume that12nX X X is a stationary 1-order Markov chain on A with transitionmatrix ()(),,k k a b a b Απ∈Π=and stationary distribution ()()a a Αµ∈. Here we are interested in the statistical distribution of the distance D between two successive occurrences of ωand more precisely in the probability()()3( occurs at 3i + 3d and there is no occurrence of at 3i + 3k, 0<k<d| occurs at 3i),D f d P D d P ωωω===where d > 0,ωoccurs at 3i means it ends at 3i . If ,3n di >− then ()0D f d =. If 13n d i≤≤−,then the event {}E= occurs at 3i + 3d ω can be decomposed into several disjoint events.{}{}()122E = occurs at 3i + 3d and there is no occurrence of at 3i + 3k, 0<k<d E = occurs at 3i + 3d and there are some occurrence of at 3i + 3k, 0<k<d E h =there is no occurrence of at 3i +3k, 0<k<h<d an ωωωωω{}d occurs at 3i + 3d and 3i+3h ωWe can get()()1| occurs at 3i .D f d PE ω=Moreover,2E can be decomposedas ()1221.d h E E h −==∪So that ()112121d h E E E E E h −===∪∪, therefore,()()()()()()()()()()()()()()()()()112133213313132131111111332332131131111| occurs at 3i | occurs at 3i | occurs at 3i ,,,,,D d h d h d d D h d d d h D h f d P E P E P E h f h f h ωωωµωµωπωωπωωµωµωµωπωωπωωµω−=−−−−=−−−−===−=−⎡⎤=−⎢⎥⎣⎦∑∑∑where 13n di≤≤−, elsewise ()0D f d =.4 Codon count distributionWe know that a codon 123ωωωω=be a word of length 3 on a finite alphabet A and ()i i X X ∈= is a random sequence on Α.This section is devoted to consider the statisticaldistribution of the count ()3N ωof the codon ωin the sequence 12n X X X . In the followingsubsection we will give an asymptotic distribution: the Gaussian regime.We assume that ()i i X X ∈= is a stationary first-order Markov chain on Α with transition probabilities (),k k a b π∈Ζand stationary distribution ()a a µ∈Α.()/331n i i N Y ω==∑, where(){}ends at position 3i in X i i Y Y I ωω==,13n i ≤≤. So we have the expectation and variance of()3N ωbelow()(){}()/3331ends at position 3i in X 3n i nE N P ωωµω===∑ ()()()()()()()()()22222333339n Var N E N E N E N ωωωωµω⎡⎤=−=−⎣⎦ ()()(){}()()()()()()()()()()2/3/3/322311133133233131212323133233212323132 ends at position 3i and 3j 32,,,32,,,33n n n i i i j n i i i j ni j j i ni j i E N E Y E Y YY nP nn n i ωµωωµωµωπωωπωωπωωµωµωπωωπωωπω==≤≠≤≤<≤−−≤<≤−⎛⎞⎛⎞⎜⎟==+⎜⎟⎜⎟⎝⎠⎜⎟⎝⎠=+=+⎛⎞=+−⋅⎜⎟⎝⎠∑∑∑∑∑()1311n i ω−=∑So we can get()()()()()()()()()()()()()231123322333212323131311232231313112,,,3392,339ni i n i i n n n Var N i n n n i ωµωµωπωωπωωπωωµωµωµωπωωµωµω−−=−−=⎛⎞⎜⎟⎝⎠⎛⎞=⎜⎟⎝⎠=+−⋅−+−⋅−∑∑ˆN N ωω− and theasymptotic variance.From section 3, we can get the maximum likelihood estimator ()3ˆN ωof ()()3E N ω: ()()()()212323322ˆN N N N ωωωωωω=Indeed, ()3N ωis a natural estimator of ()()212323,N ωωπωω, and ()()()3212323,N N ωωωπωω−is approximately a martingale as it is shown below.We introduce the martingale ()()/311|n n i i i i M Y E Y −==−∑F with ()132i i X X σ+= F . It is easyto verify that ()311|nnn EMM−=F . Moreover, we have()()(){}{}(){}()1112311223112323| ends at position 3i | ends at position 3i-1 and occurs at 3i | ends at position 3i-1 ends at 3i | ends at position 3i-1,i i i i i E Y P P E I I I ωωωωωωωωωωπωω−−−−====F F F F(){}()()()/3/311232311212323| ends at position 3i-1,,n n ii i i E Y I N ωωπωωωωπωω−====∑∑FTherefore,()()())3212323,M N N ωωωπωω=−.The next step we will establish theasymptotic normality of n M of martingales.Theorem 1: Let (),1,,n i i n ξ= be a triangular array of d-dimensional random vectors such that2,2n iE ξ<∞ , and Vbe a positive d d × matrix.Put (),,1,,,n i n n i σξξ= F ,(),,1|n i n i E ξ−F denotes the conditional expectation vector of ,n i ξand(),,1|n i n i Cov ξ−F denotes the conditional covariance matrix of ,n i ξ. If as n →∞,()(),,11|0nPn in i i i E ξ−=⎯⎯→∑F ,()(),,11|nn in i i ii Cov V ξ−=→∑F ,()(),,110,0nPn in i i iii P εξε−=∀>>⎯⎯→∑F ,Then()1,0,ni n iDN V ξ=⎯⎯→∑. ( see Ref. [9])Theorem 2: Let()()()[]212323323,1,V µωωπωωπωω=−, we have(0,)DMN V⎯⎯→as n→∞.Proof: This is an application of T heorem 1 for one-dimensional random variables()(),1|n i i i iY E Yξ−=−F()iholds for(),1|0n i iEξ−=F. We have to check that()/3,11|nn i iiVar Vξ−=→∑F as n→∞. Since i Y is a 0-1 random variable, we have()()(){}()()[]1112323323221||ends at position 3i-1,1,|i ii i i iY YIVar E E Yωωπωωπωω−−−−=−⎡⎤⎣⎦=F F FWe thus obtain()()()(){}()()[]()()()[]()()()[]/3,11/311/311/3123233231212323323212323323|||ends at position 3i-1,1,,1,,1,111,n1nn i iini i iini iiniVar VarVar Y E YVar YINnnnV asnξωωπωωπωωωωπωωπωωµωωπωωπωω−=−=−==−−−−=→∞=====→∑∑∑∑FFFThe third equation comes from ()()1|i i iE E Y EY−=F and the convergence follows from thelarge numbers:()()212212/N nωωµωω→. Finally,,2n iξ≤, so that ()2,0,4/,0n in Pεεξε∀>∀>>=, establishing condition()iii.()()())3212323,(0,)DN N N Vωωωπωω−⎯⎯→,as n→∞.Theorem 3: As n→∞()()()()()()221232333220,DN NN NNωωωωωσωω−⎞⎯⎯→⎟⎠, where()()()()()()[]()()[]2332232322212222µωσωµωµωωµωµωωµω=−−,()()()()()()321232*********,,,µωµωωπωωµωπωωπωω==,()()()32322323,µωωµωπωω=.Proof: Consider()()())3212323ˆ,nT N Nωωωπωω=−, where()()()32332322ˆ,NNωωπωωω=. We have()()()()()()()()()321232321232332332322323ˆ,,,,nnT N NNN NNωωωπωωωωπωωπωωωωωωωπωω=−−⎤⎡⎤⎦⎣⎦=−⎤⎦′=WherenM′ is the martingale()()()/3232311nn i i iiM Y E Yωωωω−=⎡⎤′=−⎣⎦∑F. We know()110,DMN V⎯⎯→ from Theorem 2, where()()()11212323323,1,V Vµωωπωωπωω⎡⎤==−⎣⎦.Similarly,()220,DnM N V′⎯⎯→, where()()()()()22/32312232332311lim,1,ni iniV Var Ynωωµωπωωπωω−→∞=⎡⎤=⎣⎦−=∑F()()[]()()()()()()(){}()()()()()(){}()()()()()(){}/3/312323111/31232311123231131111n ni j j ji jni i i i i iii i i j j jni jCov M EM M EMn n n nE M Mn nnE Y E Y Y E Yi inE Y E Y Y E YnE Y E Y Y E Ynωωωωωωωωωωωωωωωωωω−−==−−=−−≤≠≤′′′−−′=⎡⎤⎡⎤=−−∑∑⎣⎦⎣⎦⎡⎤⎡⎤=−−+∑⎣⎦⎣⎦⎡⎤⎡⎤−−⎣⎦⎣⎦∑F FF FF F()()()()()(){}/312323111ni i i i i iiE Y E Y Y E Ynωωωωωω−−=⎡⎤⎡⎤=−−∑⎣⎦⎣⎦F FThe last equation comes from()()()()()()()()()()()()()()()() ()()()()12323123231231232323220,i i i j j ji j i j j j i i i ji j i jE Y E Y Y E YE Y Y E Y E Y E Y E Y E Y YE Y Y EY E Yi jωωωωωωωωωωωωωωωωωωωωωωωω−−−−−−=−−+=−=∀≠⎡⎤⎡⎤⎣⎦⎣⎦⎡⎤⎡⎤⎡⎤⎡⎤⎣⎦⎣⎦⎣⎦⎣⎦⎡⎤⎡⎤⎣⎦⎣⎦F FF FAnd because of()()()23i i i iY Y Y Yωωωω=≡, we have()()()()()(){}()(){}()()()/31232311/3221112/311/3111112211211as n.ni i i i i iini i i i i iini i iini i iiCov E Y E Y Y E YE Y Y E Y E YnE Y E YnVar Y E YnV V V Vωωωωωω−−=−−=−=−=′⎡⎤⎡⎤−−⎣⎦⎣⎦⎡⎤=−−⎣⎦⎡⎤=−⎣⎦=−→===→∞∑∑∑F FF FFFSo we have11122122;n DnV VMNV VM⎛⎞⎛⎞⎞⎛⎞⎯⎯→⎜⎟⎜⎟⎟⎜⎟⎜⎟⎝⎠⎠⎝⎠⎝⎠′.Note that the Law of Large Numbers guarantees that almost surely()()()()2122221222NNωωωµωωµω→as n→∞.From above we can deduce thatnT converges in distribution to ()()230,Nσωwith()()()()()()()()()()()()212212323323222lim31limlim,2222221221221122122222,1,1Var T nnNMnM Mn nn NV V Vσωωωωµωωµωωµωµωµωωµωωπωωπωωµω=→∞=→∞′=→∞=+−−⎛⎞⎫⎪⎬⎪⎪⎩⎭⎛⎞⎜⎟⎝⎠⎛⎡⎤⎜⎣⎦⎝=−()()()()()()()()()()()()()()()()()()21222212323323223222123232232221222323222,1,1,µωωµωµωωπωωπωωµωµωµωµωωπωωµωµωµωµωωµωµωωµω−−⎞⎟⎜⎟⎠⎡⎤⎡⎤⎣⎦⎣⎦⎡⎤⎡⎤⎣⎦⎣⎦⎡⎤⎡⎤⎣⎦⎣⎦=−=−=−−Therefore, we have proved()()()()()()221232333220,DN NN NNωωωωωσωω−⎞⎯⎯→⎟⎠. So the proof of theorem 3 is finished.References[1] Senoussi,R.(1990) Statistique asymptotique Presque-sure de modeles statitiques conves, Ann.Inst.Henri Poincare(26):19-44.[2] Robin S, Daudin JJ (1999). Exact Distribution of word occurrences in a random sequence of letters,J. Appl. Prob(36):179-193.[3] Blom,G. and Thorburn,D.(1982). How many random digits are required until given sequences are obtained?, J.Appl. Probab(19):518-531.[4] Robin S, Daudin J.J (2001). Exact distribution of the distances between any ccurrences of a set of words. Ann. Inst. Statist. Math, 53,(4): 895-905.[5] Stefanov,V.(2003). The intersite distances between pattern occurrences in strings generated by general discrete and continuous time models: an algorithmic approach, J.Appl. Probab.: 40.[6] Gentleman, J. and Mullin, R. (1989). The distribution of the frequency of occurrence of nucleotide subsequences, based on their overlap capability, Biometrics (45): 35-52.[7] Prum,B.,Rodolphe,F., and Turckheim, et.(1995). Finding words with unexpected frequencies in dexoxyribonucleic acid sequences, J.Roy. Statist. Soc. Ser. B, (57):205-220.[8] Schbath, S. ( 1995). Compound Poisson approximation of word counts in DNA sequences, ESAIM Probab. Statist.,(1):1-16.[9] Dacunha-Castelle, D. and Dufflo, M. (1983). Probabilites at statistiques 2. Problems a temps mobile, Masson.。