BAB II
HIDDEN MARKOV MODEL
2.1. Pendahuluan Proses Stokastik dapat dipandang sebagai suatu barisan peubah acak { X t , t ∈ T } dengan T adalah parameter indeks dan X t menyatakan keadaan pada saat t. Himpunan dari semua nilai state X t yang mungkin dinamakan ruang keadaan dan biasa ditulis Ω X t . Pada umumnya T adalah bilangan bulat, di mana pada nilai negatif adlah pengamatan yang lalu. Di sini akan diambil T adalah himpunan bilangan bulat non negatif, karena kasus data yang diamati bukan berdasarkan waktu tapi lokasi. Selanjutnya, proses stokastik { X t , t ∈ T } disebut suatu proses Markov jika untuk t ∈ T berlaku
P( X t +1 = j X 0 = i0 , X 1 = i1 ,K, X t −1 = it −1 ) = P( X t +1 = j X t −1 = it −1 )
(2.1)
, dengan i0 , i1 , K , it −1 , i, j ∈ Ω X t dan P( X t +1 = j X t = i ) menyatakan peluang bersyarat suatu keadaan X t +1 jika diberikan keadaan X t , artinya peluang suatu proses berada di keadaan j pada waktu t+1 hanya bergantung pada proses yang berada di keadaan i pada waktu t. Proses Markov secara umum sudah dipaparkan dengan jelas oleh Karlin dan Taylor (1975), Gillespie (1992), dan Rogers dan William (1994). Rantai Markov adalah suatu proses Markov dengan ruang keadaan yang bernilai diskrit. Rantai markov dengan parameter diskrit adalah proses markov yang memiliki ruang keadaan diskrit dengan parameter T = (1, 2, K , n) . Teori lebih dalam tentang rantai Markov bisa dipelajari pada Ross (1996) Dalam praktek, realisasi dari rantai Markov umumnya tidak dapat diobservasi secara langsung, karena itu disebut model Markov tersembunyi. State pada model ini terbagi menjadi 2, yaitu observasi state dan state tersembunyi, sehingga model ini disebut juga model bivariat. Contoh dalam bioinformatika, jika diberikan barisan DNA, ACCAATGGATAGAC..., maka barisan DNA tersebut merupakan observasi state, dan hidden statenya bisa bermacam-macam, tergantung dari informasi apa yang ingin diperoleh. Seperti pada kasus Hongki (2004), unsur penjajaran DNA bisa dimasukkan dengan menganggap cocok (match), sisipan (insertion), dan hapusan (deletion) sebagai keadaan yang tersembunyi dengan maksud apakah barisan DNA
11
yang terobservasi tersebut berasal dari keadaan yang cocok, sisipan, atau hapusan. Contoh lain dalam bioinformatika yang bisa dijadikan sebagai keadaan yang tersembunyi adalah {intron, ekson, UTR} dan {coding region, non-coding region} untuk MMT orde-1. Untuk analisis struktur orde-2, yang tersembunyi adalah {CpGisland, non-island} dan {α-helix, β-sheet, loop}.
2.2. Definisi MMT
Model MMT merupakan 5-tuple (5 pasangan di mana masing-masing anggota bisa berupa himpunan atau ukuran) sebagai berikut: 1) Banyaknya elemen state yang tidak terobservasi (tersembunyi) pada model, N. Dianggap X t sebagai keadaan yang terjadi pada saat t. Sebagai contoh, perhatikan kasus di atas. Di sini N = 3, karena terdapat keadaan cocok, sisipan, dan hapusan. Untuk kasus ini bisa ditulis
X t = i, i = 1 (cocok), 2 (sisipan), dan 3 (hapusan). N
dapat berupa konstanta ataupun variabel acak. Hal ini sangat bergantung pada si pengamat. Dalam pembahasan di sini N adalah suatu konstanta. 2) Matrikas peluang transisi A = {aij } dimana aij adalah elemen dari A yang merupakan peluang bersyarat dari keadaan pada saat t+1, jika diketahui keadaan X pada saat t, atau aij = P( X t +1 = j X t = i ) ;1 ≤ i, j ≤ N . Karena itu A berukuran NxN. Perlu diingat bahwa aij ≥ 0 untuk setiap 1 ≤ i, j ≤ N dan
N
∑a j =1
ij
= 1 untuk
setiap 1 ≤ i ≤ N , artinya jumlah elemen setiap baris adalah 1. Untuk kasus barisan DNA dengan keadaan seperti pada 1), akan ada matriks 1 ⎛ a11
⎜ transisi A ukuran 3x3, A = 2 ⎜ a21 3 ⎜⎝ a31
a12 a22 a32
a13 ⎞ ⎟ a23 ⎟ . a33 ⎟⎠
3) Banyaknya elemen keadaan yang terobservasi, M. M umumnya tetap, ditentukan si pengamat, tetapi M juga bisa dimisalkan variabel acak. Misalkan variabel acak untuk keadaan ini adalah K, k =1,2,...,M. Untuk kasus barisan DNA, banyaknya elemen hanya 4, masing-masing A, C, G, dan T, sehingga k = 1,2,3,4. Yang menjadi permasalahan di sini adalah realisasi dari 4 elemen tersebut sepanjang T
12
4) Distribusi peluang observasi pada saat t berada pada state i (disebut juga peluang emisi), dan dinotasikan B ={bi(k)}, dimana bik = bi (k ) = P(Ot = k X t = i ) ,1 ≤ i ≤ N ,1 ≤ k ≤ M , dan K adalah observasi pada
waktu ke-t berada (bernilai) k. jadi B berukuran NxM, dan seperti di matriks transisi A, jumlah elemen setiap baris adalah 1. Secara teoritis, realisasi Ot biasa ditulis Ot = vk, 1 ≤ vk ≤ M . Hal ini dikarenakan untuk masing-masing k, vk mempunyai nilai realisasi yang berbeda. 5) State awal π 0 = {π 0 (i )} dimana π 0 (i ) = P( X 0 = i ) ,1 ≤ i ≤ N . Untuk kasus di atas,
π 0 (1) = match, π 0 (2) = insertion, dan π 0 (3) = deletion . Istilah tuple di atas berkaitan dengan himpunan dan ukuran. Di sini himpunannya diwakili oleh variabel acak. Dari definisi di atas, cukup jelas bahwa dari nilai 5-tuple (N, M, A, B, dan π0 ), terdapat 3 komponen adalah ukuran (probabilitas), yaitu A, B, dan π0 . Akibatnya MMT lebih dikenal dengan notasi λ = ( A, B, π 0 ) dengan A berukuran NxN dan B berukuran MxN. Berikut adalah ilustrasi gambar untuk model MMT O1
O2
O3
…
OT
X1
X2
X3
…
XT
Gambar 1 Ilustrasi model Markov tersembunyi
2.3. Masalah Utama dalam MMT
Secara garis besar ada 3 masalah utama yang harus diselesaikan ketika seorang peneliti menerapkan MMT di dunia nyata. Tetapi sebelumnya, perlu dipahami dengan baik 2 istilah, yaitu: 1. Orde dari suatu rantai Markov, misal rantai Markov orde-1. Di sini, pada waktu sembarang t, peluang kejadian pada satu waktu berikutnya hanya bergantung pada waktu sekarang. Lihat penjelasan untuk
13
persamaan (2.1). Sedangkan dalam rantai Markov orde-2, keadaan yang tersembunyi memiliki syarat di dua waktu sebelumnya. 2. Independensi/Kebebasan Distribusi peluang dari suatu barisan observasi O = {O1 , O2 ,..., OT } pada suatu waktu t hanya bergantung pada barisan keadaan X = { X 1 , X 2 ,..., X T } pada waktu yang sama. Dari asumsi ini, dapat didefinisikan bahwa: T
P(O X , λ A, B ,π 0 ) = ∏ P(Ot X t , λ )
(2.2)
t =1
Bukti: P(O X , λ ) = P ( O1 = o1 ,K OT = oT | X 1 = x1 ,K , X T = xT , λ )
= P ( O1 = o1 ,| X 1 = x1 ,K , X T = xT , λ ) P ( O2 = o2 ,K OT = oT | X 1 = x1 ,K , X T = xT , λ ) = P ( O1 = o1 ,| X 1 = x1 ,K , X T = xT , λ ) P ( O2 = o2 | X 1 = x1 ,K , X T = xT , λ ) P ( O3 = o3 ,K OT = oT | X 1 = x1 ,K , X T = xT , λ )
=L
= P ( O1 = o1 ,| X 1 = x1 ,K , X T = xT , λ )K P ( OT = oT ,| X 1 = x1 ,K , X T = xT , λ ) = P ( O1 = o1 ,| X 1 = x1 , λ )K P ( OT = oT ,| X T = xT , λ ) T
= ∏ P(Ot X t , λ ) (terbukti) t =1
Persamaan (2.2) sangat berguna dalam MMT, karena semua algoritma yang digunakan untuk membangun MMT memakai sifat ini. Hal ini dikarenakan di dalam persamaan tersebut kita menghitung peluang barisan observasi diketahui model, dan karena dalam MMT melibatkan keadaan tersembunyi, maka penentuan peluangnya juga bersyarat terhadap keadaan tersembunyi tersebut. Selanjutnya, berikut ini adalah tiga masalah utama yang harus diselesaikan dalam membangun MMT, yaitu: 1. Diberikan barisan observasi O = o1o2 K oT dan model λ = ( A, B, π 0 ) . Bagaimana menghitung P(O | λ ) , yaitu peluang observasi jika diberikan model. Pertanyaan ini sangat berguna jika kita bermaksud untuk memilih salah satu yang terbaik yang sesuai dengan barisan observasinya di antara beberapa model yang ada.
14
2. Diberikan barisan observasi O = o1o2 K oT dan model λ = ( A, B, π 0 ) . Bagaimana menentukan barisan state X = x1 x2 K xT yang optimal. Kesulitannya adalah bagaimana kriteria barisan state dikatakan barisan yang optimal. 3. Bagaimana menyesuaikan model parameter λ = ( A, B, π 0 ) untuk memaksimalkan
P(O | λ ) . Dengan kata lain, bagaimana mengestimasi parameter λ = ( A, B, π 0 ) . Dengan melihat kasus pada barisan DNA, misal terdapat barisan protein ACAATG, maka masalah yang pertama adalah menghitung peluang terjadinya sususan barisan DNA terobservasi, bersyarat
λ
ditetapkan/diketahui atau
P( ACAATG λ ) . Masalah yang kedua adalah mencari barisan keadaan tersembunyi yang optimal yang elemen barisannya terdiri atas cocok, sisipan, dan hapusan dari barisan observasi yang diberikan. Masalah pertama dapat diselesaikan dengan menggunakan algoritma maju mundur dan masalah yang kedua dapat diselesaikan dengan menggunakan algoritma Viterbi. Sedangkan untuk masalah yang terakhir, algoritma Baum-Welch dapat digunakan untuk mencari solusinya. 2.4. Algoritma Maju Mundur
Algoritma ini adalah proses iterasi yang didasarkan pada perhitungan peluang bersyarat melalui sifat-sifat pada peluang. Dengan menggunakan definisi peluang bersyarat, sudah bisa dihitung P(O | λ ) , tetapi operasi perhitungan yang dibutuhkan akan bertambah banyak, naik secara eksponensial, seiring dengan bertambah panjangnya barisan observasinya. Untuk menghitung P(O | λ ) diperlukan 2T .N T kali operasi, dengan N T adalah kemungkinan keadaan tersembunyi yang terjadi, jika barisan observasi adalah sepanjang T dan state sebanyak N. Ambil kasus N=3 dan M=2, dimana keadaan tersembunyinya adalah cocok (1), sisipan (2), dan hapusan (3) dan barisan DNA berukuran 2, yaitu observasi AG. Untuk kemudahan penulisan, misalkan A = 1 dan G = 3 dan ambil, ⎛ 0.7 0.1 0.2 ⎞ ⎛ 0.6 0.2 0.1 0.1 ⎞ ⎛ 0.5 ⎞ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ A0 = ⎜ 0.6 0.3 0.1 ⎟ , B0 = ⎜ 0.2 0.4 0.3 0.1 ⎟ , dan π 0 = ⎜ 0.3 ⎟ ⎜ 0.4 0.2 0.4 ⎟ ⎜ 0.1 0.1 0.5 0.3 ⎟ ⎜ 0.2 ⎟ ⎝ ⎠ ⎝ ⎠ ⎝ ⎠
15
Berdasarkan barisan observasi dan keadaan tersembunyi di atas, terlihat bahwa terdapat 9 kemungkinan yang terjadi untuk menghitung P(O1 = 1, O2 = 3 λ0 ) : 1) P ( O1 = 1, O2 = 3 dan X 1 = 1, X 2 = 1 λ ) = P ( O2 = 3, X 2 = 1 O1 = 1, X 1 = 1, λ ) .P ( O1 = 1, X 1 = 1, λ ) = ⎡⎣ P (O2 = 3 X 2 = 1, λ ). P( X 2 = 1 X 1 = 1, λ ) ⎤⎦ .P ( O1 = 1 X 1 = 1, λ ) .P( X 1 = 1, λ ) = (0.1) . (0.7) . (0.6) . (0.5) = 0.021 2) P ( O1 = 1, O2 = 3 dan X 1 = 1, X 2 = 2 λ ) = P ( O2 = 3, X 2 = 2 O1 = 1, X 1 = 1, λ ) .P ( O1 = 1, X 1 = 1, λ ) = ⎡⎣ P (O2 = 3 X 2 = 2, λ ). P( X 2 = 2 X 1 = 1, λ ) ⎤⎦ .P ( O1 = 1 X 1 = 1, λ ) .P( X 1 = 1, λ ) = (0.3) . (0.1) . (0.6) . (0.5) = 0.009 3) P ( O1 = 1, O2 = 3 dan X 1 = 1, X 2 = 3 λ ) = (0.5) . (0.2) . (0.6) . (0.5) = 0.030 4) P ( O1 = 1, O2 = 3 dan X 1 = 2, X 2 = 1 λ ) = (0.1) . (0.6) . (0.2) . (0.3) = 0.0036 5) P ( O1 = 1, O2 = 3 dan X 1 = 2, X 2 = 2 λ ) = (0.3) . (0.3) . (0.2) . (0.3) = 0.0054 6) P ( O1 = 1, O2 = 3 dan X 1 = 2, X 2 = 3 λ ) = (0.5) . (0.1) . (0.2) . (0.3) = 0.003 7) P ( O1 = 1, O2 = 3 dan X 1 = 3, X 2 = 1 λ ) = (0.1) . (0.4) . (0.1) . (0.2) = 0.0008 8) P ( O1 = 1, O2 = 3 dan X 1 = 3, X 2 = 2 λ ) = (0.3) . (0.2) . (0.1) . (0.2) = 0.0012 9) P ( O1 = 1, O2 = 3 dan X 1 = 3, X 2 = 3 λ ) = (0.5) . (0.4) . (0.1) . (0.2) = 0.004 Berdasarkan 9 nilai di atas, bisa dihitung saat ini dan besok sesuai bersyarat λ atau P(O1 = 1, O2 = 3 λ ) dengan cara menjumlahkan 9 nilai di atas,
P(O1 = 1, O2 = 3 λ ) =
∑ P (O
1
= 1, O2 = 3 dan X 1 = i, X 2 = j λ ) ,1 ≤ i, j ≤ 3
i, j
= 0.021 + 0.009 + 0.030 + 0.0036 + 0.0054 + 0.003 + 0.0008 + 0.0012 + 0.004 = 0.078
16
Perhatikan pada 1), 4), 7), 8), dan 9), terlihat bahwa peluangnya sangat kecil. Hal ini wajar karena peluang emisinya sebesar 0.1, lebih kecil dibanding yang lain. Kesimpulannya, jika nilai elemen matriks emisi sangat kecil, maka peluang yang dihasilkan juga kecil. Hal ini juga berlaku untuk matriks transisi, seperti yang ditunjukkan oleh 2) dan 6), hanya tidak sekecil 7) karena di 2) dan 6) peluang pada matriks emisi dan peluang awalnya masih cukup besar. Pada dasarnya kita bisa memperluas perumusan ini untuk barisan observasi T dengan mengambil jumlah dari peluang masing-masing kemungkinan yang terjadi. Untuk lebih jelasnya lihat persamaan (2.3) di bawah ini: P(O λ ) =
∑
P(O, X λ )
(2.3)
X = ( X1 ,..., X N )
Pada kasus di atas bisa terlihat bahwa kita butuh setidaknya 2T .N T = 2.2.32 = 36 kali operasi. Hal ini bisa dilihat dari ke-9 nilai, untuk setiap kejadian yang mungkin membutuhkan 4 kali operasi. Terlihat metode di atas tidak efisien, jadi diperlukan suatu metode baru yang dapat mereduksi operasi perhitungan P(O | λ ) . Kemudian diperkenalkan algoritma maju dan mundur. Algoritma ini menyimpan nilai yang terhitung pada iterasi yang sebelumya, sehingga paling banyak butuh N 2T operasi. Perbedaannya mungkin tidak akan terlalu berbeda secara signifikan dengan metode sebelumnya untuk kasus di atas, karena panjang barisan observasinya hanya 2. Akan tetapi algoritma ini akan sangat berguna dan efisien ketika panjang barisan observasinya cukup besar. Umumnya, T jauh lebih besar dari N, jadi pengurangan ini sangat membantu pemakai MMT. Di bawah ini dapat dipelajari algoritma maju. Definisikan α t (i ) sebagai variabel maju dimana
α t (i ) = P ( O1O2 K Ot , X t = i λ )
(2.4)
dengan kata lain, α t (i ) menyatakan total peluang observasi berakhir di state X i pada saat t dimana t = 1, 2,..., T jika diberikan barisan observasi O1O2 ...Ot . Secara umum, algoritma maju ini memuat 3 kegiatan (Rabiner 1989):
17
1. Inisialisasi
α1 ( i ) = π 0 ( i ) bi ( O1 ) untuk
i = 1,K , N
(2.5)
2. Induksi ⎧⎪ N ⎫⎪ α t +1 ( j ) = ⎨∑ α t (i )aij ⎬b j (Ot +1 ) untuk j = 1,K , N dan t = 1,K , T − 1 ⎪⎩i =1 ⎪⎭
(2.6)
3. Terminasi N
P(O | λ ) = ∑ α T (i )
(2.7)
i =1
Pada tahap inisialisasi (t=1), peluang berada di state ke i tidak lain adalah peluang awalnya, π 0 (i ) , sehingga berdasarkan (2.4) dengan t=1 diperoleh (2.5),
α1 (i ) = P ( O1 , X 1 = xi λ ) = P ( O1 | X 1 = xi , λ ) P ( X 1 = xi , λ ) = π 0 (i ) P ( O1 | X 1 = xi , λ ) , karena π 0 (i ) = P ( O1 | X 1 = xi , λ ) = π 0 (i )
bi ( O1 )
untuk i = 1,2, K , N Pada tahap induksi, akan dihitung nilai α pada saat t > 1 , pembuktiannya dilakukan dengan memanfatkan 2 asumsi MMT,
α t +1 ( j ) = P ( O1 , O2 ,..., Ot , Ot +1 , X t +1 = j | λ ) = P ( O1 , O2 ,..., Ot , Ot +1 | X t +1 = j , λ ) P( X t +1 = j | λ ) = P ( O1 , O2 ,..., Ot | X t +1 = j , λ ) P(Ot +1 | X t +1 = j , λ ) P( X t +1 = j | λ ) = P ( O1 , O2 ,..., Ot , X t +1 = j | λ ) P(Ot +1 | X t +1 = j , λ ) = P ( O1 , O2 ,..., Ot , X t +1 = j | λ ) b j (Ot +1 ) ⎡N ⎤ = ⎢ ∑ P ( O1 , O2 ,..., Ot , X t = i, X t +1 = j λ ) ⎥ b j (Ot +1 ) ⎣ i =1 ⎦ ⎡N ⎤ = ⎢ ∑ P ( O1 , O2 ,..., Ot , X t = i λ ) P( X t +1 = j | O1 , O2 ,..., Ot , X t = i, λ) ⎥ b j (Ot +1 ) ⎣ i =1 ⎦ ⎡N ⎤ = ⎢ ∑ P ( O1 , O2 ,..., Ot , X t = i λ ) P( X t +1 = j | X t = i, λ) ⎥ b j (Ot +1 ) ⎣ i =1 ⎦ ⎡N ⎤ = ⎢ ∑ α t (i )aij ⎥ b j (Ot +1 ) ⎣ i =1 ⎦
dengan j = 1,K , N dan t = 1,K , T − 1
18
Sebagai contoh, misalkan hanya ada 2 observasi (t=2), maka untuk menghitung α 2 ( j ) , bisa dilihat langkah-langkahnya sebagai berikut:
α 2 ( j ) = P ( O1O2 , X 2 = j λ ) = P ( O1O2 , X 1 = 1, X 2 = j λ ) + P ( O1O2 , X 1 = 2, X 2 = j λ ) + K + P ( O1O2 , X 1 = N , X 2 = j λ ) N
= ∑ P ( O1O2 , X 1 = i, X 2 = j λ ) i =1 N
= ∑ P ( O1O2 | X 1 = i, X 2 = j , λ ) P ( X 1 = i, X 2 = j , λ ) i =1 N
= ∑ P ( O1O2 | X 1 = i, X 2 = j , λ ) P ( X 2 = j | X 1 = i, λ ) P ( X 1 = i, λ ) i =1 N
= ∑ π 0 ( i ) aij P ( O1O2 | X 1 = i, X 2 = j , λ ) i =1 N
= ∑ π 0 ( i ) aij bi ( O1 ) b j ( O2 ) i =1
N
= ∑ α1 (i )aij b j (O2 ) i =1
⎧⎪ N ⎫⎪ = b j (O2 )⎨∑ α1 (i )aij ⎬ ⎪⎩i =1 ⎪⎭ Jadi dengan induktif dapat ditulis secara umum sebagai berikut
α t +1 ( j ) = ⎧⎨∑ α t (i )aij ⎫⎬b j (Ot +1 ) untuk N
⎩i =1
⎭
j = 1,K , N
dan t = 1,K , T − 1
(2.8)
Maka N
P(O | λ ) = ∑ α T (i ) i =1
(2.9)
Maksud dari (2.9) adalah menjumlahkan semua peluang gabungan dari observasi dan keadaan tersembunyi bersyarat model sehingga diperoleh peluang marjinal dari observasi tersebut. Sebagai ilustrasi dari algoritma maju, perhatikan gambar berikut,
19
α t −1 (1)
a1
a1,1b1 ( Ot )
,2
b1
(O t
)
α t +1 (1) α t (1)
α t +1 ( 2 )
α t −1 ( 2 )
αt ( 2 ) α t −1 ( 3)
α t +1 ( 3)
α t ( 3)
Gambar 2. Ilustrasi Algoritma maju
Sebagai contoh, ambil kasus pada metode sebelumnya (9 nilai peluang) yang hasilnya selanjutnya dibandingkan. Akan dihitung P(O1 = 1, O2 = 3 λ ) dari suatu barisan observasi AG. Perhatikan kembali matriks matriks transisi, matriks emisi, dan nilai awal pada hal. 15, Pada tahap inisialisasi, tentu ada 3 buah α , masing-masing untuk cocok, sisipan, dan hapusan,
α1 (1) = π (1) b1 (1) = (0.5) . (0.6) = 0.3 α1 ( 2 ) = π ( 2 ) b2 (1) = (0.3) . (0.2) = 0.06 α1 ( 3) = π ( 3) b3 (1) = (0.2) . (0.1) = 0.02 Selanjutnya, pada tahap induksi akan dihitung α 2 ( j ) dengan 1 ≤ j ≤ 3 . ⎧
3
⎫
α 2 (1) = ⎨∑ α1 ( i ) ai1 ⎬ b1 ( O2 = 3) ⎩ i =1
⎭
= {(0.3) . (0.7) + (0.06) . (0.6) + (0.02) . (0.4)} . (0.1) = 0.0254 ⎧
3
⎫
α 2 ( 2 ) = ⎨∑ α1 ( i ) ai 2 ⎬ b2 ( O2 = 3) ⎩ i =1
⎭
= {(0.3) . (0.1) + (0.06) . (0.3) + (0.02) . (0.2)} . (0.3) = 0.0156
20
⎧
3
⎫
α 2 ( 3) = ⎨∑ α1 ( i ) ai 3 ⎬ b3 ( O2 = 3) ⎩ i =1
⎭
= {(0.3) . (0.2) + (0.06) . (0.1) + (0.02) . (0.4)} . (0.5) = 0.037 Terakhir adalah tahap terminasi dengan cara menghitung P(O1 = 1, O2 = 3 λ ) dengan menjumlahkan α 2 ( j ) dengan 1 ≤ j ≤ 3 , P(O1 = 1, O2 = 3 λ )
=
3
∑α (i ) i =1
2
= 0.0254 + 0.0156 + 0.037 = 0.078 Berdasarkan contoh ini, terlihat dua metode memiliki hasil yang sama, namun banyaknya operasi perhitungan tidak sama. Terlihat bahwa algoritma maju lebih efisien daripada metode definisi peluang bersyarat. Selanjutnya, akan dibahas metode mundur beserta contohnya di kasus yang sama. Hal ini dilakukan untuk mengetahui apakah metode maju dan mundur menghasilkan peluang yang sama. Di bawah ini adalah metode mundur. Algoritma ini langkahnya serupa dengan algoritma maju. Di sini inisialisasi didasarkan pada seluruh observasi yang ada, jadi algoritma ini mengganti O1O2 K Ot pada persamaan (2.4) dengan Ot +1Ot + 2 K OT . Di sini pemakai berpikir seperti akan menaksir ke depan. Pandang variabel mundur,
β t (i ) = P(Ot +1Ot + 2 K OT | X t = i , λ )
(2.10)
Algoritma mundur selengkapnya dapat dilihat di bawah ini: 1. Inisialisasi
β T (i ) = 1 , untuk i = 1,K , N
2. Induksi
(2.11)
N
βt ( i ) = ∑ b j (Ot +1 )βt +1 ( j )aij , untuk t = T − 1,T − 2,K ,1 , i = 1,K , N
(2.12)
j =1
3. Terminasi
21
N
P ( O | λ ) = ∑ bi (1)π 0 (i ) β1 ( i )
(2.13)
i =1
Pada tahap inisialisasi, diambil β T (i ) = 1 karena i adalah state final, dan bernilai nol untuk i yang lain. Untuk tahap induksi, pembuktiannya mirip dengan algoritma maju. Untuk lebih jelasnya, langkah-langkah pembuktiannya bisa dilihat di bawah ini:
βt ( i )
= P ( Ot +1 , Ot + 2 ,..., OT X t = i, λ ) N
= ∑ P ( Ot +1 , Ot + 2 ,..., OT , X t +1 = j X t = i, λ ) j =1 N
= ∑ P ( Ot +1 Ot + 2 ,..., OT , X t +1 = j , X t = i, λ )P ( Ot + 2 ,..., OT , X t +1 = j X t = i, λ ) j =1 N
= ∑ P ( Ot +1 X t +1 = j , λ )P ( Ot + 2 ,..., OT X t +1 = j , X t = i, λ ) P ( X t +1 = j X t = i, λ ) j =1 N
= ∑ P ( Ot +1 X t +1 = j , λ )P ( Ot + 2 ,..., OT X t +1 = j , λ ) P ( X t +1 = j X t = i, λ ) j =1 N
= ∑ P ( Ot +1 X t +1 = j , λ )β t +1 ( j ) P ( X t +1 = j X t = i, λ ) j =1 N
= ∑ b j (Ot +1 )β t +1 ( j )aij j =1
untuk t = T − 1,T − 2,K ,1 , i = 1,K , N Sebagai contoh, kita ambil kasus yang sama dengan kasus pada algoritma maju untuk menghitung P(O1 = 1, O2 = 3 λ ) . Pada tahap inisialisasi, β 2 ( i ) = 1 untuk i = 1, 2, 3. Selanjutnya, di tahap induksi akan dihitung β1 ( i ) dengan i = 1, 2, 3,
β1 (1) =
3
∑a j =1
1j
b j (3) β 2 ( j )
= (0.7) . (0.1) . 1 + (0.1) . (0.3) . 1 + (0.2) . (0.5) . 1 = 0.2
β1 ( 2 ) =
3
∑a j =1
2j
b j (3) β 2 ( j )
= (0.6) . (0.1) . 1 + (0.3) . (0.3) . 1 + (0.1) . (0.5) . 1 = 0.2
22
β1 ( 3) =
3
∑a j =1
3j
b j (3) β 2 ( j )
= (0.4) . (0.1) . 1 + (0.2) . (0.3) . 1 + (0.4) . (0.5) . 1 = 0.3 Untuk memulai induksi, dibandingkan dengan algoritma maju di mana proporsi α didominasi oleh keadaan cocok dan hapusan, di sini proporsi β untuk masing-masing keadaan hampir sama. Jadi pada λ yang sama tidak ada keterkaitan nilai awal α1 dan β1 . Pada tahap terminasi, akan dihitung P(O1 = 1, O2 = 3 λ ) dengan P(O1 = 1, O2 = 3 λ ) =
3
∑ b (1)π β (i) i =1
i
i
1
= (0.6) . (0.5) . (0.2) + (0.2) . (0.3) . (0.2) + (0.1) . (0.2) . (0.3) = 0.06 + 0.012 + 0.006 = 0.078 Berdasarkan contoh kasus di atas, terlihat bahwa perhitungan melalui algoritma maju dan mundur akan menghasilkan nilai peluang yang sama. Tapi kesimpulan ini tidak selalu benar. Jika elemen-elemen peluang transisi dan emisi berupa bilangan rasional, seperti contoh di atas, maka algoritma maju dan mundur akan memberikan hasil yang sama. Sebaliknya, jika matriks tersebut mengandung satu atau lebih bilangan irasional, maka metode maju dan mundur akan memberikan hasil yang berbeda. Dengan kata lain, jika peluang transisi dan emisi ditaksir dari frekuensi relatif, metode maju dan mundur yang dihasilkan sering kali berbeda. Variabel mundur bisa juga dikombinasikan dengan variabel maju untuk menghitung P ( O λ ) seperti yang tampak di bawah ini: P ( O1 ,..., OT , X t = i λ ) = P ( O1 ,..., OT X t = i, λ ) P( X t = i, λ )
= P ( O1 ,..., Ot X t = i, λ ) P ( Ot +1 ,..., OT X t = i, λ ) P( X t = i, λ ) = P ( O1 ,..., Ot , X t = i λ ) P ( Ot +1 ,..., OT X t = i, λ )
23
= α t (i ) β t (i ) N
P ( O λ ) = ∑ α t (i ) βt (i )
(2.14) (2.15)
i =1
2.5. Algoritma Viterbi
Setelah mempelajari algoritma maju dan mundur untuk menghitung P ( O λ ) , akan dilanjutkan untuk menyelesaikan masalah kedua dalam membangun MMT, yaitu jika diberikan barisan observasi O = o1o2 K oT dan model λ = ( A, B ,π ) , bagaimana menentukan barisan hidden state X = x1 x2 K xT yang optimal. Salah satu metode yang sering digunakan adalah dengan algoritma Viterbi. Definisikan,
δt (i ) =
max
X1 , X 2 ,K, X t −1
P ( O1O2 K Ot , X 1 X 2 L X t −1 , X t = i λ )
(2.16)
Berikut ini adalah algoritma lengkap untuk Viterbi: 1. Inisialisasi
δ1 ( i ) = π 0 (i )bi ( O1 ) , 1 ≤ i ≤ N ψ 1 (i ) = 0
(2.18)
δ t ( j ) = b j (Ot ) max {aijδ t −1 (i )} , 2 ≤ t ≤ T dan
(2.19)
ψ t ( j ) = arg max {aijδ t −1 ( i )} , 2 ≤ t ≤ T dan 1 ≤ j ≤ N
(2.20)
(2.17)
2. Induksi 1≤i ≤ N
1≤i ≤ N
3. Terminasi P* = max {δ T (i )}
(2.21)
qT* = arg max{δT (i )}
(2.22)
1≤ i ≤ N
1≤i ≤ N
4. Penjajagan Mundur
( )
q*t = ψ t +1 q*t+1 , t = T − 1,T − 2,L ,1
(2.23)
Pada tahap inisialisasi (t = 1),
24
δ1 ( i ) = P ( X 1 = i, O1 ) = P ( O1 | X 1 = i ) P ( X 1 = i ) = bi ( O1 ) π 0 (i ) Pada tahap induksi,
δt ( j) =
max
X1 , X 2 ,K, X t −1
P ( O1 K Ot , X 1 L X t −1 , X t = j λ )
= max
X1 , X 2 ,K, X t −1
=
{P ( O
t
max
}
O1 K Ot −1 , X 1 L X t −1 , X t = j , λ ) P ( O1 K Ot −1 , X 1 L X t −1 , X t = j , λ )
X1 , X 2 ,K, X t −1
{P ( O
}
X t = j , λ ) P ( O1 K Ot −1 , X 1 L X t −1 , X t = j , λ )
t
= P ( Ot X t = j , λ ) = b j (Ot )
max
max
max { P ( O1 K Ot −1 , X 1 L X t − 2 , X t −1 = i, X t = j , λ )}
X1 , X 2 ,K, X t −2 1≤i ≤ N
{
{
= b j (Ot ) max P ( X t = j X t −1 = i ) 1≤i ≤ N
}
max P ( X t = j X t −1 = i ) P ( O1 K Ot −1 , X 1 L X t − 2 , X t −1 = i, λ )
X1 , X 2 ,K, X t −2 1≤i ≤ N
{
max
X1 , X 2 ,K, X t −2
}
P ( O1 K Ot −1 , X 1 L X t − 2 , X t −1 = i, λ )
}
= b j (Ot ) max P ( X t = j X t −1 = i ) δ t −1 (i ) 1≤i ≤ N
= b j (Ot ) max {aijδ t −1 (i )} 1≤i ≤ N
Sebagai contoh, ambil kasus untuk barisan protein dengan observasi AG sama seperti pada kasus-kasus sebelumnya. Berdasarkan algoritma Viterbi, diperoleh: 1. Inisialisasi
δ1 (1) = π 0 (1)b1 (1) = (0.5) . (0.6) = 0.3 δ1 ( 2 ) = π 0 (2)b2 (1) = (0.3) . (0.2) = 0.06 δ1 ( 3) = π 0 (3)b3 (1) = (0.2) . (0.1) = 0.02 2. Induksi
δ 2 (1) = b1 (3) max {ai1δ1 (i )} 1≤i ≤3
= (0.1) . max{(0.7) . (0.3),(0.6) . (0.06),(0.4) . (0.02)} = (0.1) . (0.7) . (0.3) = 0.021
ψ 2 (1) = 1
25
δ 2 (2) = b2 (3) max {ai 2δ1 (i )} 1≤ i ≤ 3
= (0.3) . max{(0.1) . (0.3),(0.3) . (0.06),(0.2) . (0.02)} = (0.3) . (0.1) . (0.3) = 0.009
ψ 2 (2) = 1 δ 2 (3) = b3 (3) max {ai 3δ1 (i )} 1≤ i ≤ 3
= (0.5) . max{(0.2) . (0.3),(0.1) . (0.06),(0.4) . (0.02)} = (0.5) . (0.2) . (0.3) = 0.03
ψ 2 (3) = 1 3. Terminasi P* = max {δ 2 ( i )} 1≤ i ≤3
= max {δ 2 (1) , δ 2 ( 2 ) , δ 2 ( 3)} = max{0.021, 0.009, 0.03} = 0.03 q2* = arg max {δ 2 ( i )} 1≤i ≤3
=3 4. Penjajagan Mundur q1* = ψ 2 ( q2* ) = ψ 2 ( 3) = 1
Berdasarkan hasil algoritma Viterbi di atas, dapat disimpulkan bahwa barisan hidden state yang optimal adalah x1 = cocok dan x3 = hapusan . 2.6. Algoritma Baum-Welch
Untuk mengatasi masalah ketiga dalam MMT, biasanya digunakan algoritma Baum-Welch. Algoritma Baum-Welch digunakan untuk mengestimasi 3 parameter MMT
sehingga terbentuk model baru λˆ ( Aˆ , Bˆ , πˆ ) di mana P(O λˆ) ≥ P(O λ ) .
Prosedur ini diulang terus menerus sampai proses konvergen. Permasalahan yang muncul adalah bagaimana menentukan nilai awal π 0 agar barisan taksiran tersebut konvergen dengan cepat, sehingga iterasi yang dilakukan efisien. Pemilihan nilai awal
26
ini perlu dilakukan dengan hati-hati, karena hasil dari taksiran parameter MMT sangat dipengaruhi oleh nilai awal ini. Algoritma Baum-Welch merupakan kasus khusus dari algoritma EM (Ekspektasi Maksimum). Algoritma EM merupakan algoritma yang dipakai untuk mempelajari model-model probabilistik dalam suatu kasus yang melibatkan keadaan yang tersembunyi (Rabiner 1989). Lebih khusus, algoritma ini merupakan prosedur iteratif untuk menghitung penaksir Maksimum-likelihood, L(λ ) , yang diterapkan pada data tidak lengkap. Algoritma EM memiliki 2 prosedur berurutan yang dilakukan secara berulang secara terpadu. 2 prosedur tersebut adalah: 1. E-step : menentukan ekspektasi bersyarat E [log P(O, X | λ )] 2. M-step : memaksimumkan ekspektasi bersyarat tersebut terhadap λ Dalam kasus MMT, data lengkap adalah informasi mengenai keadaan Xt sama baiknya dengan observasi Ot untuk setiap t. Jadi data tidak lengkapnya adalah hanya observasi Ot itu sendiri. Fungsi likelihood dari data tak lengkap adalah P(O | λ ) dan fungsi likelihood dari data lengkap adalah P(O, X | λ ) . Pada E-step, data yang “hilang/tersembunyi” ditaksir melalui data yang terobservasi dan estimasi parameter model. Peluangnya di-log-kan untuk mempermudah perhitungan (G. F. Larrahondo, S. Bridges, Eric A. H., 2005). Untuk peluang yang sangat kecil, kekonvergenan akan lebih cepat tercapai dan membutuhkan batas ambang yang juga sangat kecil. Untuk menghindari hal-hal ini, maka peluang perlu di-log-kan. Dengan kata lain, dalam algoritma EM, yang akan diestimasi adalah E [log P(O, X | λ )]
(2.24)
Berdasarkan definisi ekspektasi, persamaan (2.24) dapat ditulis sebagai berikut,
∑ P ( X | O, λ ) log ( P ( O, x | λ ) ) 0
(2.25)
x
dengan λ0 merupakan taksiran awal dari parameter λ . Berdasarkan asumsi kebebasan bersyarat, maka diperoleh T
P ( O, X | λ ) = π 0 ∏ axt −1 xt bxt ( Ot ) t =1
(2.26)
dengan T menyatakan banyaknya observasi (panjang observasi). Kemudian dengan mensubsitusikan persamaan (2.26) ke persamaan (2.25), diperoleh
27
⎛ T ⎞ Q ( λ , λ0 ) = ∑ P ( O, x | λ0 ) log π 0 + ∑ P ( O, x | λ0 ) ⎜ ∑ log axt −1 xt ⎟ x x ⎝ t =1 ⎠ T ⎛ ⎞ + ∑ P ( O, x | λ0 ) ⎜ ∑ log bxt ( Ot ) ⎟ x ⎝ t =1 ⎠ Selanjutnya perhatikan bagian pertama pada persamaan (2.27), N
∑ P ( O, x | λ ) log π = ∑ P ( O, X 0
0
i =1
x
0
= i | λ0 ) log π 0 (i )
dengan menggunakan multiplier Langrange dengan syarat bahwa
(2.27)
(2.28) N
∑π i =1
0
(i ) = 1 ,
diperoleh ⎛ N ⎛ N ⎞⎞ P O , X = i | λ log π ( i ) + η π 0 (i) − 1⎟ ⎟ = 0 ( ) ∑ ⎜∑ 0 0 0 ⎜ ∂π 0 (i ) ⎝ i =1 ⎝ i =1 ⎠⎠ ∂
maka P ( O, X 0 = i | λ0 ) +η = 0 π 0 (i)
π 0 (i) = −
P ( O, X 0 = i | λ0 )
η
N
∑π i =1
N
(i ) = −∑
0
1= −
(2.29)
P ( O, X 0 = i | λ0 )
η
i =1
1
η
N
∑ P ( O, X i =1
0
= i | λ0 )
N
η = −∑ P ( O, X 0 = i | λ0 )
(2.30)
i =1
sedangkan, N ∑ P(O, X 0 i =1
N
= i | λ0 ) = ∑ P(O, X 0 = i | λ0 )P( X 0 = i, λ0 ) i =1
= P(O | λ0 )
dari persamaan (2.29) dan (2.30) maka diperoleh,
π 0 (i ) = −
P ( O, X 0 = i | λ0 )
=
P ( O, X 0 = i | λ0 )
P ( O | λ0 ) Berikutnya adalah untuk bagian kedua dari persamaan (2.27) adalah
η
⎛T ∑ P (O, x | λ 0 )⎜ ∑ log a x x t −1 t x ⎝ t =1
⎞ ⎟ ⎠
(2.31)
(2.32) N
Dengan menggunakan multiplier Langrange dengan syarat ∑ aij = 1 , j =1
28
∂ ∂aij
⎛ ⎛ T P O , x | λ ( ) ⎜⎜ ∑ 0 ⎜ ∑ log a xt −1 xt x ⎝ t =1 ⎝
⎛ N ⎞⎞ ⎞ + η ⎜ ∑ aij − 1⎟ ⎟⎟ = 0 ⎟ ⎠ ⎝ j =1 ⎠⎠
∂ ∂aij
⎛ N N ⎛ T ⎜⎜ ∑∑ P ( O, X t −1 = i, X t = j | λ0 ) ⎜ ∑ log axt −1 xt ⎝ t =1 ⎝ i =1 j =1
∂ ∂aij
⎛ N N T ⎛ N ⎞⎞ ⎜⎜ ∑∑∑ log aij P ( O, X t −1 = i, X t = j | λ0 ) + η ⎜ ∑ aij − 1⎟ ⎟⎟ = 0 ⎝ j =1 ⎠⎠ ⎝ i =1 j =1 t =1
maka
T
∑ P ( O, X t =1
t −1
= i, X t = j | λ0 ) aij
T
aij = −
∑ P ( O, X t =1
t −1
⎛ N ⎞⎞ ⎞ η + ⎜ ∑ aij − 1⎟ ⎟⎟ = 0 ⎟ ⎠ ⎝ j =1 ⎠⎠
+η
= 0 atau
= i, X t = j | λ0 )
, terlihat bahwa aij masih bergatung
η
pada η , untuk itu dicari nilai η dengan menjumlahkan aij untuk semua j, diperoleh T
N
∑a j =1
ij
N
= −∑
∑ P ( O, X t =1
η
j =1
N
= i, X t = j | λ0 )
t −1
atau
T
η = −∑∑ P ( O, X t −1 = i, X t = j | λ0 ) , mengingat j =1 t =1
P ( O, X t −1 = i, X t = j | λ0 ) = P ( O, X t −1 = i | X t = j , λ0 ) P ( X t = j , λ0 ) , maka T
N
η = −∑∑ P ( O, X t −1 = i | X t = j, λ0 ) P ( X t = j, λ0 ) atau t =1 j =1 T
η = −∑ P ( O, X t −1 = i | λ0 ) t =1
Selanjutnya η disubtitusi ke aij di atas, sehingga diperoleh ∑ P (O, X t −1 T
aij =
t =1
= i, X t = j | λ0 )
∑ P (O, X t −1 T
t =1
= i | λ0 )
(2.33)
Selanjutnya, perhatikan bagian ketiga pada persamaan (2.27), yaitu ∑ P (O, x | λ0 )⎛⎜ ∑ log bxt (ot )⎞⎟ T
x ⎝ t =1 ⎠ Dengan cara yang sama dengan sebelumnya, menggunakan multiplier Langrange
dengan syarat ∑ bi ( j ) =1, didapat K
j =1
29
⎛ N ⎛ K ⎞⎞ ⎛ T ⎞ P O , X i | λ log b o η = + ( ) ( ) ⎜⎜ ∑ ⎜ ∑ bi ( j ) − 1⎟ ⎟⎟ = 0 t i t ⎟ 0 ⎜∑ ∂bi ( j ) ⎝ i =1 ⎝ t =1 ⎠ ⎝ j =1 ⎠⎠ ∂
atau ⎛ N N T ⎛ K ⎞⎞ ⎜⎜ ∑∑∑ log bi ( j ) P ( O, X t = i | λ0 ) δ j + η ⎜ ∑ bi ( j ) − 1⎟ ⎟⎟ = 0 ∂bi ( j ) ⎝ j =1 i =1 t =1 ⎝ j =1 ⎠⎠ ∂
dengan ⎧1 , ot = j ⎩0 , lainnya
δj =⎨ Setelah diturunkan diperoleh, N
T
∑∑ P ( O, X i =1 t =1
t
= i | λ0 ) δ j
bi ( j ) N
bi ( j ) = −
T
∑∑ P ( O, X i =1 t =1
K
∑ b ( j ) = −∑ j =1
i
N
= i | λ0 ) δ j
(2.34)
T
∑∑ P ( O, X i =1 t =1
j =1
K
atau
η N
K
t
+η = 0
t
= i | λ0 ) δ j
η
T
N
T
η = −∑∑∑ P ( O, X t = i | λ0 ) δ j = ∑∑ P ( O, X t = i | λ0 ) j =1 i =1 t =1
i =1 t =1
dengan mensubstitusikan η pada persamaan (2.34), diperoleh, ∑ ∑ P (O, X t N T
bi ( j ) =
i =1 t =1 N T
= i | λ0 )δ j
∑ ∑ P(O, X t
i =1 t =1
= i | λ0 )
∑ P(O, X t T
=
t =1 T
= i | λ0 )δ j
∑ P (O, X t
t =1
= i | λ0 )
(2.35)
Algoritma Baum-Welch juga dikenal sebagai algoritma Maju-Mundur dengan variabel maju dan mundur adalah sebagai berikut: ⎪⎧ α t (i ) = P ( O1O2 KOt , X t = i λ ) ⎨ ⎪⎩ β t ( i ) = P ( Ot +1Ot + 2 KOT , X t = i | λ )
(2.36)
Kemudian didefinisikan variabel baru
ξt (i, j ) = P ( X t = i, X t +1 = j O, λ )
(2.37)
Dengan menggunakan definisi peluang bersyarat dan aturan Bayes, maka persamaan (2.37) berubah menjadi
30
ξt (i, j ) = P ( X t = i, X t +1 = j O, λ ) =
P ( X t = i, X t +1 = j , O λ ) P (O λ )
= P ( O1O2 K Ot , X t = i λ ) P ( X t +1 = j X t = i ) P ( Ot +1 X t +1 = j ) P ( Ot + 2 K OT , X t +1 = j | λ ) P (O λ )
=
α t (i)aij b j (Ot +1 ) βt +1 ( j ) P (O λ )
dengan P ( O λ ) =
(2.38) N N −1
N
∑ α (i) β (i) = ∑∑ α (i)a b (O i =1
t
t
t
i =1 j =1
ij
t +1
j
) β t +1 ( j ) .
Dari ξt (i, j ) yang diketahui, bisa dihitung peluang berada di state i pada saat t, γ t (i ) , dengan menjumlahkan atas j, N
∑ ξ (i, j )
γ t (i) =
j =1
(2.39)
t
Dihitung ekspektasi banyaknya transisi dari state i ke state j, sebut saja dengan
τ (i, j ) , dengan menjumlahkan atas t, T
τ (i, j ) = ∑ ξt (i, j )
(2.40)
t =1
Kemudian dihitung ekspektasi banyaknya yang berpindah dari state i, tulis dengan
τ (i ) , dengan menjumlahkan atas t dari 1 sampai T, T
T
N
τ (i ) = ∑ γ t (i) = ∑∑ ξt (i, j ) t =1
(2.41)
t =1 j =1
Dengan menghubungkan antara persamaan (2.31) dengan (2.39), maka πˆi bisa berbentuk:
πˆ0 (i) =
P ( O, X 0 = i | λ0 ) P ( O | λ0 )
= P ( X 0 = i O , λ0 ) N −1
=
∑ P( X j =1
0
= i , X 1 = j O, λ )
31
N −1
=
∑ ξ (i, j ) j =1
0
= γ 0 (i )
(2.42)
Kemudian Aˆ = aˆij bisa juga diubah dalam bentuk lain dengan menghubungkan persamaan (2.33), (2.40), dan (2.41): T
aˆij =
∑ P ( O, X t =1
t −1
T
∑ P ( O, X t =1
T −1
=
∑ P ( O, X t =0
t
T −1
∑ P ( O, X T −1
∑ P( X t =0
t
T −1
∑P(X T −1
∑ P( X t =0
t
t
= i | λ)
t
= i | O, λ0 )P(O λ )
= i, X t +1 = j | O, λ )
T −1
∑ P( X t =0
=
= i |λ)
= i, X t +1 = j | O, λ ) P(O λ )
t =0
=
t −1
= i, X t +1 = j | λ )
t =0
=
= i, X t = j | λ )
t
= i | O, λ0 )
τ (i, j ) τ (i)
(2.43)
Terakhir, Bˆ = bˆi ( j ) bisa ditulis dalam bentuk yang lebih sederhana dengan melihat persamaan (2.35), (2.39), dan (2.41): T
bˆi ( j ) =
∑ P ( O, X t =1 T
∑ P ( O, X t =1
T
=
t
∑ P( X t =1 T
t
∑P(X t =1
= i | λ )δ j t
= i | λ)
= i | O, λ ) δ j t
= i | O, λ )
32
T
∑
=
t =1,Ot = j T
γ t ( j) (2.44)
∑ γ ( j) t =1
t
Sebagai contoh, kita tinjau kembali kasus barisan protein AG yang pernah dibahas pada bagian-bagian sebelumnya,
ξˆ1 (1,1) =
α1 (1)a11b1 (3) β 2 (1) (0.3) . (0.7) . (0.1) . (1) = = 0.269 0.078 P(O λ )
ξˆ1 (1, 2) =
α1 (1)a12b2 (3) β 2 (2) (0.3) . (0.1) . (0.3) . (1) = = 0.115 P(O λ ) 0.078
ξˆ1 (1,3) =
α1 (1)a13b3 (3) β 2 (3) (0.3) . (0.2) . (0.5) . (1) = = 0.385 P(O λ ) 0.078
ξˆ1 (2,1) =
α1 (2)a21b1 (3) β 2 (1) (0.06) . (0.6) . (0.1) . (1) = = 0.046 P(O λ ) 0.078
sama seperti di atas, dan berturut-turut diperoleh ξˆ1 (2, 2) = 0.069 , ξˆ1 (2,3) = 0.039 ,
ξˆ1 (3,1) = 0.011 , ξˆ1 (3, 2) = 0.015 , dan ξˆ1 (3,3) = 0.051 . γˆ1 (1) = γˆ1 (2) = γˆ1 (3) =
3
∑ ξˆ (1, j) = 0.269 + 0.115 + 0.385 = 0.769 j =0
1
3
∑ ξˆ (2, j) = 0.046 + 0.069 + 0.039 = 0.154 j =0
1
3
∑ ξˆ (3, j) = 0.011 + 0.015 + 0.051 = 0.077 j =0
1
2 −1
Karena hanya ada 2 observasi, maka τˆ(1) = ∑ γˆt (1) = γˆ1 (1) = 0.769, τˆ(2) = γˆ1 (2) = t =1
0.154, dan τˆ(3) = γˆ1 (3) = 0.077. Nilai τˆ(i, j ) = ξˆ1 (i, j ) untuk 1 ≤ i, j ≤ 3 . Dari hasil perhitungan di atas, dapat kita estimasi 3 parameter MMT, ⎛ γˆ1 (1) ⎞ ⎛ 0.769 ⎞ πˆ0 (i) = ⎜⎜ γˆ1 (2) ⎟⎟ = ⎜⎜ 0.154 ⎟⎟ ⎜ γˆ (3) ⎟ ⎜ 0.077 ⎟ ⎝ 1 ⎠ ⎝ ⎠
33
⎛ τˆ(1,1) τˆ(1, 2) τˆ(1,3) ⎞ ⎛ 0.269 ⎜ ˆ τˆ(1) τˆ(1) ⎟⎟ ⎜ 0.769 ⎜ τ (1) ⎜ ⎜ τˆ(2,1) τˆ(2, 2) τˆ(2, 3) ⎟ ⎜ 0.046 aˆij = ⎜ ⎟= τˆ(2) τˆ(2) ⎟ ⎜ 0.154 ⎜ τˆ(2) ⎜ ⎜ τˆ(3,1) τˆ(3, 2) τˆ(3, 3) ⎟ ⎜ 0.011 ⎟ ⎜ ⎜ τˆ(3) τˆ(3) τˆ(3) ⎠ ⎝ 0.077 ⎝ ⎛ 1 ⎜ ∑ γˆ1 (1) ⎜ t =1,Ot =1 ⎜ 1 ⎜ ∑ γˆ1 (1) ⎜ t =1 ⎜ 1 ⎜ ∑ γˆ1 (2) t =1,O =1 bˆij = ⎜ 1 t ⎜ ⎜ ∑ γˆ1 (2) ⎜ t =1 ⎜ 1 ⎜ ∑ γˆ1 (3) ⎜ t =1,1Ot =1 ⎜ γˆ1 (3) ⎜ ∑ = t 1 ⎝
1
∑
t =1,Ot = 2 1
γˆ1 (1)
∑ γˆ (1) 1
t =1 1
∑
γˆ1 (2)
t =1,Ot = 2 1
∑ γˆ (2) t =1
1
1
∑
t =1,Ot = 2 1
γˆ1 (3)
∑ γˆ (3) t =1
1
1
∑
0.385 ⎞ 0.769 ⎟ ⎛ 0.35 0.15 0.5 ⎞ ⎟ 0.039 ⎟ ⎜ ⎟ = ⎜ 0.3 0.45 0.25 ⎟ ⎟ 0.154 ⎟ ⎜ 0.14 0.2 0.66 ⎟⎠ 0.051 ⎟ ⎝ ⎟ 0.077 ⎠
0.115 0.769 0.069 0.154 0.015 0.077
γˆ1 (1)
1
∑
⎞
γˆ1 (1) ⎟
⎟ ⎟ γˆ1 (1) γˆ1 (1) ⎟ ∑ ∑ t =1 t =1 ⎟ 1 1 ⎟ γˆ1 (2) ∑ γˆ1 (2) ⎟ ⎛ 1 0 0 0 ⎞ ∑ t =1,Ot = 3 t =1,Ot = 4 ⎟ = ⎜1 0 0 0 ⎟ 1 1 ⎟ ⎟ ⎜ γˆ1 (2) γˆ1 (2) ⎟ ⎜⎝ 1 0 0 0 ⎟⎠ ∑ ∑ t =1 t =1 ⎟ 1 1 ⎟ γˆ1 (3) γˆ1 (3) ⎟ ∑ ∑ t =1,Ot =3 t =1,Ot = 4 ⎟ 1 1 ⎟ γˆ1 (3) γˆ1 (3) ⎟ ∑ ∑ t =1 t =1 ⎠ t =1,Ot =3 1
t =1,Ot = 4 1
Dari estimator yang dihasilkan, maka diperoleh model λˆ ( Aˆ , Bˆ , πˆ 0 ) . Dengan
( )
bentuk Bˆ = bˆij seperti di atas, maka sudah bisa dipastikan bahwa P O λˆ =1 dan ini jauh lebih besar dibandingkan dengan P ( O λ ) = 0.078 . Contoh di atas merupakan contoh yang tidak bagus, karena hanya melibatkan dua obervasi, akibatnya estimator yang dihasilkan, terutama matriks peluang emisi, akan memiliki elemen matriks yang isinya 0 dan 1. Dengan demikian maka nilai peluang observasinya akan bernilai 1, sehingga tidak perlu ada iterasi lagi.
34
