BAB III METODE JUKES-CANTOR DAN METODE KIMURA
Pada Bab ini akan dibahas mengenai metode Jukes-Cantor dan metode Kimura dalam membentuk suatu evolutionary model. Pada dasarnya kedua metode menggunakan ide yang sama dalam menggambarkan substitusi pada single site, yaitu bahwa penggambaran substitusi nukleotida pada single site ini berjalan sesuai dengan suatu regular continuous-time Markov chain [5]. Perbedaan utama antara keduanya terletak pada pembentukan matriks probabilitas transisi P ( t ) dari Markov chain yang nantinya digunakan untuk menentukan jarak antara dua barisan. Berikut dijelaskan lebih lanjut mengenai kedua metode.
3.1 METODE JUKES-CANTOR
Metode Jukes-Cantor menghasilkan suatu model yang biasa disebut sebagai model satu parameter. Model ini mengasumsikan tingkat substitusi yang sama antara transisi dengan transversi. Atau dengan kata lain, tingkat substitusi suatu site berubah dari suatu state menjadi state lain dalam S = { A,C,G,T } adalah sama. Sehingga untuk sembarang state dalam S ,
probabilitas transisi antar state-state tersebut juga sama, misal sama dengan 31 Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
32
f ( t ) . Dengan demikian, dengan diasumsikan berjalan sesuai dengan suatu
regular continuous-time Markov chain, maka pada Doron-Faigenboim, Adi [2] ditentukan matriks probabilitas transisi dari Jukes-Cantor sebagai berikut A
C
G
T
1 − 3f ( t ) f (t ) f (t ) f (t ) A 1 − 3f ( t ) f (t ) f (t ) f (t ) C , P (t ) = f (t ) f (t ) 1 − 3f ( t ) f (t ) G f (t ) f (t ) 1 − 3f ( t ) T f ( t )
(3.1)
dengan f ( t ) merupakan fungsi bernilai non negatif yang diferensiabel untuk setiap t ≥ 0 . Pada subbab ini akan dijelaskan bentuk elemen-elemen dari matriks P ( t ) (3.1) yang pembahasannya diberikan sebagai berikut.
Karena Jukes-Cantor diasumsikan berjalan sesuai dengan regular Markov chain maka seluruh elemen P ( t ) dari Jukes-Cantor diferensiabel untuk setiap t ≥ 0 sehingga elemen-elemen matriks P ( t ) dari Jukes-Cantor tersebut dapat diturunkan terhadap t sebagai berikut
−3f ' ( t ) f ' ( t ) f ' (t ) f ' (t ) f ' (t ) f ' (t ) f ' ( t ) −3f ' ( t ) P '(t ) = . f ' (t ) f ' (t ) −3 ' f ( t ) f ' ( t ) −3f ' ( t ) f ' (t ) f ' (t ) f ' ( t ) Saat t = 0 ,
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
33
−3f ' ( 0 ) f ' ( 0 ) f ' (0) f ' ( 0) −3f ' ( 0 ) f ' (0) f ' (0) f ' ( 0) P '(0) = . f ' (0) f ' (0) −3 ' f ( 0 ) f ' ( 0 ) f ' (0) f ' (0) −3f ' ( 0 ) f ' ( 0 ) Misal didefinisikan suatu matrix of instantaneous change sebagai berikut −3α α Q = P ' (0) = α α
α −3α α α
α α −3α α
α α , α −3α
dimana α = f ' ( 0 ) merupakan suatu konstanta tingkat substitusi yang bernilai positif. Berdasarkan lemma pada subbab 2.4.1.2, untuk t ≥ 0 berlaku P ' (t ) = P (t ) P ' (0) .
Maka
f ' (t ) f ' (t ) −3f ' ( t ) f ' ( t ) −3f ' ( t ) f ' (t ) f ' (t ) f ' (t ) f ' (t ) f ' ( t ) −3 ' f ( t ) f ' ( t ) f ' (t ) f ' (t ) −3f ' ( t ) f ' ( t ) 1 − 3f ( t ) f (t ) f (t ) f ( t ) −3α f (t ) 1 − 3f ( t ) f (t ) f ( t ) α ⋅ = f (t ) 1 − 3f ( t ) f (t ) f (t ) α 1 − 3f ( t ) α f (t ) f (t ) f ( t )
α −3α α α
α α −3α α
α α . α −3α
Pandang elemen pada baris ke satu kolom ke dua dari perkalian matriks tersebut, f ' ( t ) = (1 − 3f ( t ) ) α + f ( t )( −3α ) + f ( t ) ⋅ α + f ( t ) ⋅ α
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
34
f ' ( t ) = α − 3αf ( t ) − 3αf ( t ) + αf ( t ) + αf ( t )
d (f (t )) dt
d (f (t ))
α − 4αf ( t )
= α − 4αf ( t )
= dt .
Integralkan kedua ruas maka
d (f (t )) ∫ α − 4αf ( t ) = ∫ dt ln ( α − 4αf ( t ) ) −4α
= t + c1
ln ( α − 4αf ( t ) ) = −4αt − 4αc1 α − 4αf ( t ) = e −4 αt −4 αc1
= e −4 αt e −4αc1 = ce −4αt −4αf ( t ) = ce −4αt − α f (t ) =
1 c −4 αt e . − 4 4α
(3.2)
Anggap nilai awal f ( 0 ) = 0 . Sehingga untuk t = 0 , 1 c −4α ⋅0 − =0, e 4 4α
f (0) =
diperoleh nilai c = α . Substitusikan nilai c ke (3.2) diperoleh f (t ) =
1 1 −4αt − e 4 4
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
(3.3)
35
Dengan mensubstitusikan (3.3) ke (3.1) maka diperoleh elemen-elemen matriks probabilitas transisi P ( t ) dari Jukes-Cantor sebagai berikut 1 3 −4αt 4 + 4 e , i = j pij ( t ) = 1 − 1 e −4 αt , i ≠ j 4 4
(3.4)
untuk setiap i , j ∈ S , dan t ≥ 0 . Berikutnya akan ditunjukkan dua kondisi yang harus dipenuhi oleh Jukes-Cantor dalam menggambarkan substitusi nukleotida. 1)
Stationary probability distribution φ yang bersifat unik Untuk mencari stationary probability distribution φ = ( φ A , φC , φG , φT )
dari Jukes-Cantor, perhatikan sistem persamaan linier yang terbentuk dari −3α α φQ = [ φA φC φG φT ] α α
α −3α α α
α α −3α α
α α = [0 0 0 0 ] = 0 α −3α
berikut: (1)
−3φA + φC + φG + φT = 0
(2)
φA − 3φC + φG + φT = 0
(3)
φA + φC − 3φG + φT = 0
(4)
φA + φC + φG − 3φT = 0
Solusi untuk sistem persamaan linier tersebut adalah φA = φC = φG = φT .
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
36
Akan dicari φA = φC = φG = φT yang memenuhi kondisi
∑φ i∈S
i
= φA + φC + φG + φT = 1 dengan φi ≥ 0 untuk setiap i ∈ S .
Substitusikan φA = φC = φG = φT ke
∑φ i∈S
i
= φA + φC + φG + φT = 1 diperoleh
φA = φC = φG = φT =
1 . 4
Berdasarkan Definisi 1 pada subbab 2.4.1.2 maka stationary probability distribution dari Jukes-Cantor adalah
1 1 1 1 φ = , , , . 4 4 4 4
(3.5)
Ini merupakan satu-satunya stationary probability distribution dari JukesCantor sehingga Jukes-Cantor memiliki stationary probability distribution φ yang bersifat unik. 2)
Model yang bersifat time-reversible Akan dibuktikan bahwa Jukes-Cantor merupakan model yang bersifat
time-reversible. Bukti. Misalkan π = (π A ,π C ,π G ,πT ) merupakan suatu vektor probabilitas awal
1 1 1 1 dari Jukes-Cantor dengan π = φ = , , , . Berdasarkan Definisi 2 pada 4 4 4 4 subbab 2.4.1.2, maka reversed Markov chain dari Jukes-Cantor didefinisikan oleh
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
37
pij * ( t ) =
φ j p ji ( t ) φi
,
untuk setiap i , j ∈ S , dan t ≥ 0 . Karena φi = φ j =
1 untuk setiap i , j ∈ S maka 4 pij * ( t ) = p ji ( t ) ,
untuk setiap i , j ∈ S , dan t ≥ 0 . Berdasarkan elemen-elemen matriks probabilitas transisi P ( t ) dari JukesCantor (3.4), maka untuk setiap i , j ∈ S , dan t ≥ 0 diperoleh 1 3 −4αt 1 3 −4 αt 4 + 4 e , j = i 4 + 4 e , i = j pij * ( t ) = p ji ( t ) = = = pij ( t ) . 1 1 1 1 − 4 αt − 4 αt − e , j ≠i − e , i≠ j 4 4 4 4 Atau dengan kata lain P * ( t ) = P ( t ) untuk t ≥ 0 . Sehingga berdasarkan Definisi 3 pada subbab 2.4.1.2, model Jukes-Cantor merupakan model yang bersifat time-reversible. Berdasarkan 1) dan 2) maka kedua kondisi yang harus dipenuhi oleh Jukes-Cantor terpenuhi. Model Jukes-Cantor diperoleh dengan mengestimasi jarak antara dua barisan berbeda dengan panjang yang sama menggunakan maximum likelihood. Berdasarkan model likelihood (2.8) pada subbab 2.4.2, L = ln φx1 + ln φx2 + ... + ln φxn + ln px1y1 ( t ) + ln px2 y 2 ( t ) + K + ln pxn y n ( t ) .
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
38
Pandang elemen-elemen matriks probabilitas transisi P ( t ) (3.4) dan stationary probability distribution φ (3.5) dari Jukes-Cantor. Maka model likelihood untuk Jukes-Cantor dapat ditulis menjadi lnL = ln
1 1 1 + ln + ... + ln + m1 ln psama ( t ) + m2 ln pbeda ( t ) 4 4 4
= n ln
1 + m1 ln psama ( t ) + m2 ln pbeda ( t ) 4
(3.6)
dimana m1 + m2 = n , dengan m1
= banyaknya substitusi ke state yang sama,
m2
= banyaknya substitusi ke state yang berbeda,
psama ( t ) = probabilitas site berubah ke state yang sama, pbeda ( t ) = probabilitas site berubah ke state yang berbeda,
untuk t ≥ 0 . Dengan memaksimumkan model likelihood (3.6) maka akan diperoleh
d (ln L ) dt
=0
1 d n ln + m1 ln psama ( t ) + m2 ln pbeda ( t ) 4 =0 dt m1
psama ' ( t ) psama ( t )
+ m2
pbeda ' ( t ) pbeda ( t )
=0
(3.7)
Berdasarkan matriks probabilitas transisi P ( t ) dari Jukes-Cantor (3.4) maka
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
39
psama ( t ) =
d ( pii ( t ) ) 1 3 −4 αt + e = −3αe −4 αt , dan ⇒ psama ' ( t ) = 4 4 dt
(3.8)
pbeda ( t ) =
d ( pij ( t ) ) 1 1 −4αt − e ⇒ pbeda ' ( t ) = = αe −4 αt . 4 4 dt
(3.9)
Substitusikan (3.8) dan (3.9) ke (3.7), maka model likelihood (3.7) sekarang menjadi m1
−3αe −4 αt αe −4 αt + m2 =0 1 1 −4αt 1 3 −4 αt + e − e 4 4 4 4
Karena α merupakan suatu konstanta yang bernilai positif maka m1
−3e −4αt e −4αt + m2 =0. 1 1 −4αt 1 3 −4 αt + e − e 4 4 4 4
Misal z = e −4αt , maka
m1
−3 z z + m2 =0 1 1 1 3 + z − z 4 4 4 4
−3m1
z z + m2 =0 1− z 1 + 3z
−3m1 (1 − z ) z + m2 (1 + 3z ) z = 0
−3m1z + 3m1z 2 + m2 z + 3m2 z 2 = 0
( 3m1 + 3m2 ) z 2 − ( 3m1 − m2 ) z = 0 diperoleh
z=
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
3m1 − m2 . 3m1 + 3m2
(3.10)
40
Dengan mensubstitusikan kembali nilai z = e −4αt maka
e −4αt =
3m1 − m2 3m1 + 3m2
−4αt = ln
t =−
dimana r =
3m1 − m2 3m1 + 3m2 3m1 − m2 1 ln 4α 3m1 + 3m2
=−
1 3m1 − m2 + 3m2 − 3m2 ln 4α 3m1 + 3m2
=−
1 ( 3m1 + 3m2 ) − ( −m2 − 3m2 ) ln 4α 3m1 + 3m2
=−
4m2 1 ln 1 − 4α 3m1 + 3m2
=−
1 4 m2 ln 1 − 4α 3 m1 + m2
=−
1 4 ln 1 − r , 4α 3
m2 m banyaknya substitusi ke state yang berbeda = 2 = . m1 + m2 n panjang barisan
Karena pada subbab 2.4.2 t didefinisikan sebagai jarak antara dua barisan ( d ) maka diperoleh
d=−
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
1 4 ln 1 − r . 4α 3
41
Dengan pemilihan nilai α =
1 maka diperoleh model Jukes-Cantor yang 3
umum digunakan untuk menghitung jarak antara dua barisan yang berbeda sebagai berikut
3 4 d = − ln 1 − r 4 3 dimana r =
(3.11)
banyaknya substitusi ke state yang berbeda . panjang barisan
3.2 METODE KIMURA
Metode Kimura menghasilkan suatu model yang biasa disebut sebagai Kimura 2-Parameter model atau model dua parameter. Model ini merupakan generalisasi dari model Jukes-Cantor yang telah dibahas pada subbab sebelumnya. Model ini mengasumsikan tingkat substitusi yang berbeda antara transisi dengan transversi sehingga probabilitas transisi antara substitusi transisi dengan substitusi transversi tersebut juga berbeda, misal masing-masing sama dengan f ( t ) dan g (t ) . Sehingga dengan diasumsikan berjalan sesuai dengan suatu regular continuous-time Markov chain maka pada Penn, Osnat [8] ditentukan matriks probabilitas transisi dari Kimura sebagai berikut
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
42
A
C
G
T
1 − f ( t ) − 2g ( t ) A g (t ) f (t ) g (t ) g (t ) f (t ) 1 − f ( t ) − 2g ( t ) g (t ) C P (t ) = G 1 − f ( t ) − 2g ( t ) g (t ) g (t ) f (t ) 1 − f ( t ) − 2g ( t ) T g (t ) f (t ) g (t )
(3.12) dengan f ( t ) dan g (t ) merupakan fungsi-fungsi bernilai non negatif yang diferensiabel untuk setiap t ≥ 0 . Sama seperti pada pembahasan model Jukes-Cantor sebelumnya, akan dicari bentuk elemen-elemen dari matriks P ( t ) (3.12). Karena Kimura diasumsikan berjalan sesuai dengan regular Markov chain maka seluruh elemen P ( t ) dari Kimura diferensiabel untuk setiap t ≥ 0 sehingga elemen-elemen matriks P ( t ) dari Kimura tersebut dapat diturunkan terhadap t sebagai berikut
g ' (t ) f ' (t ) g ' (t ) −f ' ( t ) − 2g ' ( t ) g ' (t ) g ' (t ) f ' (t ) −f '(t ) − 2g ' ( t ) . P ' (t ) = f ' (t ) g ' (t ) g ' (t ) −f ' ( t ) − 2g ' ( t ) g ' (t ) g ' (t ) f ' (t ) −f '(t ) − 2g ' ( t ) Saat t = 0 ,
g ' (0) f ' (0) g ' (0) −f ' ( 0 ) − 2g ' ( 0 ) g ' (0) g ' (0) f ' (0) −f ' ( 0 ) − 2g ' ( 0 ) . P ' (0) = f ' (0) g ' (0) g ' (0) −f ' ( 0 ) − 2g ' ( 0 ) g ' (0) f ' (0) g ' (0) −f ' ( 0 ) − 2g ' ( 0 ) Misal didefinisikan suatu matrix of instantaneous change sebagai berikut
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
43
β α β −α − 2 β β −α − 2 β β α , Q = P ' (0) = α β −α − 2 β β α β −α − 2 β β dimana α = f ' ( 0 ) dan β = g ' ( 0 ) merupakan konstanta-konstanta tingkat substitusi yang bernilai positif, dengan α merupakan tingkat substitusi dari transisi dan β merupakan tingkat substitusi dari transversi. Berdasarkan lemma pada subbab 2.4.1.2, untuk t ≥ 0 berlaku P ' (t ) = P (t ) P ' (0) .
Maka
g ' (t ) f ' (t ) g ' (t ) −f ' ( t ) − 2g ' ( t ) g ' (t ) g ' (t ) f ' (t ) −f '(t ) − 2g ' ( t ) f ' (t ) g ' (t ) g ' (t ) −f ' ( t ) − 2g ' ( t ) g ' (t ) f ' (t ) −f '(t ) − 2g ' ( t ) g ' (t ) g (t ) f (t ) g (t ) 1 − f ( t ) − 2g ( t ) g (t ) g (t ) f (t ) 1 − f ( t ) − 2g ( t ) ⋅ = f (t ) g (t ) g (t ) 1 − f ( t ) − 2g ( t ) 1 − f ( t ) − 2g ( t ) g (t ) g (t ) f (t ) β α β −α − 2 β β −α − 2 β β α . α β −α − 2 β β α β −α − 2 β β Pandang elemen pada baris ke satu kolom ke tiga dari perkalian matriks tersebut, f ' ( t ) = (1 − f ( t ) − 2g ( t ) ) α + g ( t ) β + f ( t )( −α − 2β ) + g ( t ) β
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
44
f ' ( t ) = α − αf ( t ) − 2αg ( t ) + βg ( t ) − αf ( t ) − 2βf ( t ) + βg ( t )
d (f (t )) dt
= α − 2 (α + β ) f ( t ) − 2 ( α − β ) g ( t ) .
(3.13)
Pandang kembali elemen pada baris ke satu kolom ke dua dari perkalian matriks tersebut, g ' ( t ) = (1 − f ( t ) − 2g ( t ) ) β + g ( t )( −α − 2 β ) + f ( t ) β + g ( t ) α = β − βf ( t ) − 2βg ( t ) − αg ( t ) − 2βg ( t ) + βf ( t ) + αg ( t )
d (g (t )) dt
d (g (t ))
β − 4 βg ( t )
= β − 4 βg ( t )
= dt .
Integralkan kedua ruas maka
d ( g (t ))
∫ β − 4 βg ( t ) = ∫ dt ln ( β − 4 βg ( t ) ) −4 β
= t + c1
ln ( β − 4 βg ( t ) ) = −4 βt − 4 βc1 β − 4 βg ( t ) = e −4 βt − 4 βc1
= e −4 βt e −4 βc1 = ce −4 βt −4βg ( t ) = ce −4 βt + β
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
45
g (t ) =
1 c −4 βt e . − 4 4β
(3.14)
Anggap nilai awal g ( 0 ) = 0 . Sehingga untuk t = 0 ,
g (0) =
1 c −4 β⋅0 =0, − e 4 4β
diperoleh nilai c = β . Substitusikan nilai c ke (3.14) diperoleh g (t ) =
1 1 −4 βt − e . 4 4
(3.15)
Substitusikan (3.15) ke (3.13) diperoleh
d (f (t ))
1 1 = α − 2 ( α + β ) f ( t ) − 2 ( α − β ) − e −4 βt 4 4
dt
= α − 2 (α + β ) f (t ) − = −2 ( α + β ) f ( t ) +
d (f (t )) dt
+ f (t ) (2 (α + β )) =
α − β α − β −4 βt e + 2 2
α + β α − β −4 βt + e 2 2
α + β α − β −4 βt + e . 2 2
(3.16)
Persamaan diferensial (3.16) ini merupakan suatu persamaan diferensial linier. Untuk menyelesaikannya, kalikan kedua ruas dengan faktor integrasi
e∫
2( α + β )dt
.
d (f (t )) 2(α + β )dt α + β α − β −4 βt ∫ 2(α + β )dt + f (t ) ( 2 (α + β )) e∫ e e + = dt 2 2
d ( f ( t ) ) ∫ 2( α + β )dt 2 ( α + β )dt α + β ∫ 2( α + β )dt α − β −4 βt ∫ 2( α + β )dt e e e e = + + f (t ) ( 2 (α + β )) e∫ dt 2 2
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
46
d ( f ( t ) ) ∫ 2( α + β )dt e dt
d e ∫ + f (t )
2( α + β )dt
dt
α + β 2( α + β )dt α − β 2( α + β )dt = e∫ e −4 βt e ∫ + 2 2
d f ( t ) e ∫ dt
p( t )dt
α + β 2( α + β )dt α − β 2( α + β )dt = + e∫ e −4 βt e ∫ 2 2
d f ( t ) e ∫
p ( t )dt
= α + β e ∫ 2( α + β )dt + α − β e −4 βt e ∫ 2( α + β )dt dt 2 2
Integralkan kedua ruas maka f (t ) e∫
2( α + β )dt
f (t ) e
2( α + β )t
α + β α − β −4 βt ∫ 2( α + β )dt = ∫ + dt e e 2 2 α + β α − β −4 βt 2(α + β )t dt e e + = ∫ 2 2 α + β 2( α + β )t α − β 2( α − β )t = ∫ e + e dt 2 2 =
α−β e ( ) α+β e ( ) ⋅ + +c 2 2 (α + β ) 2 2 (α − β )
=
1 2( α + β )t 1 2( α − β )t +c e + e 4 4
f (t ) =
2 α +β t
2 α−β t
1 1 −4 βt −2 α + β t + e + ce ( ) 4 4
Anggap nilai awal f ( 0 ) = 0 . Sehingga untuk t = 0 , f (0) =
1 1 −4 β ⋅0 −2 α + β ⋅0 + ce ( ) = 0 , + e 4 4
1 diperoleh nilai c = − . Substitusikan nilai c ke (3.17) diperoleh 2
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
(3.17)
47
f (t ) =
1 1 −4 βt 1 −2( α + β )t + e − e 4 4 2
(3.18)
Dengan mensubstitusikan (3.15) dan (3.18) ke (3.12) maka diperoleh elemen-elemen matriks probabilitas transisi P ( t ) dari Kimura sebagai berikut (1) pAA ( t ) = pCC ( t ) = pGG ( t ) = pTT ( t ) =
1 1 −4 βt 1 −2( α + β )t + e + e , 4 4 2
(2) pAG ( t ) = pGA ( t ) = pCT ( t ) = pTC ( t ) =
1 1 −4 βt 1 −2( α + β )t + e − e , dan 4 4 2
(3) pAC ( t ) = pCA ( t ) = pAT ( t ) = pTT ( t ) = pCG ( t ) = pGC ( t ) = pGT ( t ) = pTG ( t )
=
1 1 −4 βt − e . 4 4
(3.19)
untuk t ≥ 0 . Berikutnya akan ditunjukkan dua kondisi yang harus dipenuhi oleh Kimura dalam menggambarkan substitusi nukleotida. 1)
Stationary probability distribution φ yang bersifat unik Untuk mencari stationary probability distribution φ = ( φ A , φC , φG , φT )
dari Kimura, perhatikan sistem persamaan linier yang terbentuk dari β α β −α − 2 β β −α − 2 β β α φQ = [ φA φC φG φT ] = [0 0 0 0 ] α β −α − 2 β β α β −α − 2 β β =0 berikut:
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
48
(1)
( −α − 2 β ) φA + βφC + αφG + βφT
(2)
βφ A + ( −α − 2 β ) φC + βφG + αφT = 0
(3)
αφA + βφC + ( −α − 2 β ) φG + βφT = 0
(4)
βφ A + αφC + βφG + ( −α − 2 β ) φT = 0
=0
Solusi untuk sistem persamaan linier tersebut adalah φA = φC = φG = φT . Akan dicari φA = φC = φG = φT yang memenuhi kondisi
∑φ i∈S
i
= φA + φC + φG + φT = 1 dengan φi ≥ 0 untuk setiap i ∈ S .
Substitusikan φA = φC = φG = φT ke
∑φ i∈S
i
= φA + φC + φG + φT = 1 diperoleh
φA = φC = φG = φT =
1 . 4
Berdasarkan Definisi 1 pada subbab 2.4.1.2 maka stationary probability distribution dari Kimura adalah
1 1 1 1 φ = , , , . 4 4 4 4
(3.20)
1 1 1 1 φ = , , , ini merupakan satu-satunya stationary probability distribution 4 4 4 4 dari Kimura sehingga Kimura memiliki stationary probability distribution φ yang bersifat unik.
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
49
2)
Model yang bersifat time-reversible Akan dibuktikan bahwa Kimura merupakan model yang bersifat time-
reversible. Bukti. Misalkan π = (π A ,π C ,π G ,πT ) merupakan suatu vektor probabilitas awal
1 1 1 1 dari Kimura dengan π = φ = , , , . Berdasarkan Definisi 2 pada 4 4 4 4 subbab 2.4.1.2, maka reversed Markov chain dari Kimura didefinisikan oleh pij * ( t ) =
φ j p ji ( t ) φi
,
untuk setiap i , j ∈ S , dan t ≥ 0 . Karena φi = φ j =
1 untuk setiap i , j ∈ S maka 4 pij * ( t ) = p ji ( t ) ,
untuk setiap i , j ∈ S , dan t ≥ 0 . Berdasarkan elemen-elemen matriks probabilitas transisi P ( t ) dari Kimura (3.22), maka untuk t ≥ 0 diperoleh (1) pAA * ( t ) = pCC * ( t ) = pGG * ( t ) = pTT * ( t ) =
1 1 −4 βt 1 −2( α + β )t + e + e , 4 4 2
(2) pAG * ( t ) = pGA * ( t ) = pCT * ( t ) = pTC * ( t ) =
1 1 −4 βt 1 −2( α + β )t + e − e , dan 4 4 2
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
50
(3) pAC * ( t ) = pCA * ( t ) = pAT * ( t ) = pTT * ( t ) = pCG * ( t ) = pGC * ( t ) = pGT * ( t ) = pTG * ( t ) =
1 1 −4 βt − e . 4 4
Atau dengan kata lain P * ( t ) = P ( t ) untuk t ≥ 0 . Sehingga berdasarkan Definisi 3 pada subbab 2.4.1.2, model Kimura merupakan model yang bersifat time-reversible. Berdasarkan 1) dan 2) maka kedua kondisi yang harus dipenuhi oleh Kimura terpenuhi. Seperti yang telah disebutkan sebelumnya, model Kimura merupakan generalisasi dari model Jukes-Cantor. Hal ini dikarenakan pada tingkat substitusi transisi yang sama dengan tingkat substitusi transversi, yaitu pada pemilihan nilai α = β , maka model Kimura akan memiliki bentuk yang sama dengan Jukes-Cantor. Berikut diberikan penjelasannya. Berdasarkan model likelihood (2.8) pada subbab 2.4.2, L = ln φx1 + ln φx2 + ... + ln φxn + ln px1y1 ( t ) + ln px2 y 2 ( t ) + K + ln pxn y n ( t ) .
Pandang elemen-elemen matriks probabilitas transisi P ( t ) (3.19) dan stationary probability distribution φ (3.20) dari Kimura. Maka model likelihood untuk Kimura dapat ditulis menjadi lnL = ln
1 1 1 + ln + ... + ln + m1 ln psama ( t ) + m2 ln ptransisi ( t ) + m3 ln ptransversi ( t ) 4 4 4
= n ln
1 + m1 ln psama ( t ) + m2 ln ptransisi ( t ) + m3 ln ptransversi ( t ) , 4
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
(3.21)
51
dimana m1 + m2 + m3 = n , dengan m1
= banyaknya substitusi ke state yang sama,
m2
= banyaknya substitusi transisi,
m3
= banyaknya substitusi transversi,
psama ( t )
= probabilitas site berubah ke state yang sama,
ptransisi ( t ) = probabilitas substitusi transisi, ptransversi ( t ) = probabilitas substitusi transversi,
untuk t ≥ 0 . Dengan memaksimumkan model likelihood (3.21) maka akan diperoleh
d (ln L ) dt
=0
1 d n ln + m1 ln psama ( t ) + m2 ln ptransisi ( t ) + m3 ln ptransversi ( t ) 4 =0 dt m1
psama ' ( t ) psama ( t )
+ m2
ptransisi ' ( t ) ptransisi ( t )
+ m3
ptransversi ' ( t ) ptransversi ( t )
=0.
(3.22)
Berdasarkan matriks probabilitas transisi P ( t ) dari Kimura (3.19) maka psama ( t ) =
1 1 −4 βt 1 −2( α + β )t −2 α + β t + e + e ⇒ psama ' ( t ) = − βe −4 βt − ( α + β ) e ( ) (3.23) 4 4 2
ptransisi ( t ) =
1 1 −4 βt 1 −2( α + β )t −2 α + β t − e + e ⇒ ptransisi ' ( t ) = − βe −4 βt + ( α + β ) e ( ) (3.24) 4 4 2
ptransversi ( t ) =
1 1 −4 βt − e ⇒ ptransversi ' ( t ) = βe −4 βt . 4 4
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
(3.25)
52
Substitusikan (3.23), (3.24), dan (3.25) ke (3.22), maka model likelihood (3.22) sekarang menjadi − βe −4 βt + ( α + β ) e ( ) − βe −4 βt − ( α + β ) e ( ) βe −4 βt + m2 + m3 =0 1 1 −4 βt 1 −2(α + β )t 1 1 −4 βt 1 −2(α + β )t 1 1 −4 βt + e + e + e − e − e 4 4 2 4 4 2 4 4 −2 α + β t
−2 α + β t
m1
(3.26) dimana α dan β merupakan konstanta-konstanta tingkat substitusi yang bernilai positif dengan α merupakan tingkat substitusi dari transisi dan β merupakan tingkat substitusi dari transversi. Pilih α = β , maka model (3.26) menjadi m1
αe −4 αt −αe −4αt − 2αe −4αt −αe −4αt + 2αe −4αt + m2 + m3 =0 1 1 −4 αt 1 −4 αt 1 1 −4 αt 1 −4 αt 1 1 −4 αt + e + e + e − e − e 4 4 2 4 4 2 4 4 −3αe −4 αt αe −4αt αe −4αt m1 + m2 + m3 =0 1 3 −4 αt 1 1 −4αt 1 1 −4αt + e − e − e 4 4 4 4 4 4 m1
−3αe −4 αt αe −4αt + ( m2 + m3 ) =0, 1 3 −4 αt 1 1 −4αt + e − e 4 4 4 4
atau m1
−3αe −4 αt αe −4 αt + m4 =0, 1 3 −4 αt 1 1 −4 αt + e − e 4 4 4 4
(3.27)
dimana m1 adalah banyaknya substitusi ke state yang sama dan m4 = m2 + m3 adalah banyaknya substitusi ke state yang berbeda. Model (3.27) ini merupakan model yang sama seperti yang diberikan oleh model
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008
53
maximum likelihood Jukes-Cantor (3.10). Sehingga dengan kata lain, untuk
α = β , Kimura akan menghasilkan model yang sama dengan model JukesCantor. Berikut diberikan model Kimura untuk menghitung jarak antara dua barisan.
(
)
1 d = − ln (1 − 2r1 − r2 ) 1 − 2r2 , 2
(3.28)
dimana r1 =
banyaknya substitusi transisi banyaknya substitusi transversi , dan r2 = . panjang barisan panjang barisan Perolehan model Kimura (3.28) tersebut tidak akan dijelaskan lebih
lanjut pada skripsi ini, namun dapat dilihat perolehan selengkapnya pada Kimura, Motoo [6]. Model Jukes-Cantor (3.11) dan model Kimura (3.28) ini selanjutnya digunakan untuk menghitung jarak tiap pasang barisan yang akan disimulasikan dalam bentuk program pada Bab IV terhadap suatu data barisan DNA.
Suatu Kajian..., Febrini Cesarina, FMIPA UI, 2008