Algoritma Expectation-Maximization (EM) merupakan metode iteratif yang digunakan untuk memperoleh estimasi maximum likelihood pada model statistik yang melibatkan variabel laten (latent variable) yang tidak dapat diamati secara langsung.
Dalam simulasi ini, terdapat dua koin, yaitu Koin A dan Koin B, dengan probabilitas kemunculan sisi Head yang belum diketahui (\(\theta_A\) dan \(\theta_B\)). Eksperimen dilakukan melalui serangkaian pelemparan koin tanpa mengetahui identitas koin yang digunakan pada setiap set percobaan.
Data observasi terdiri dari 5 set percobaan, di mana setiap set melibatkan 10 kali pelemparan koin:
# Menyiapkan data observasi
data_H <- c(5, 9, 8, 4, 7)
n_lemparan <- 10
set_names <- paste("Set", 1:5)Pada simulasi pertama, nilai awal parameter ditetapkan sebesar \(\theta_A = 0.6\) dan \(\theta_B = 0.5\).
Tahapan yang dilakukan:
# Inisialisasi
theta_A <- 0.6
theta_B <- 0.5
iterasi <- 10
log_ppt <- data.frame(Iterasi = 0:iterasi, Theta_A = NA, Theta_B = NA)
log_ppt[1, 2:3] <- c(theta_A, theta_B)
for(i in 1:iterasi) {
# E-Step: Menghitung responsibilitas (P(Z|E))
p_E_ZA <- dbinom(data_H, n_lemparan, theta_A)
p_E_ZB <- dbinom(data_H, n_lemparan, theta_B)
mu_A <- p_E_ZA / (p_E_ZA + p_E_ZB)
mu_B <- p_E_ZB / (p_E_ZA + p_E_ZB)
# M-Step: Update parameter
theta_A <- sum(data_H * mu_A) / sum(mu_A * n_lemparan)
theta_B <- sum(data_H * mu_B) / sum(mu_B * n_lemparan)
log_ppt[i+1, 2:3] <- c(round(theta_A, 3), round(theta_B, 3))
}
knitr::kable(log_ppt, caption = "Hasil Iterasi dengan Parameter Awal Standar")| Iterasi | Theta_A | Theta_B |
|---|---|---|
| 0 | 0.600 | 0.500 |
| 1 | 0.713 | 0.581 |
| 2 | 0.745 | 0.569 |
| 3 | 0.768 | 0.550 |
| 4 | 0.783 | 0.535 |
| 5 | 0.791 | 0.526 |
| 6 | 0.795 | 0.522 |
| 7 | 0.796 | 0.521 |
| 8 | 0.796 | 0.520 |
| 9 | 0.797 | 0.520 |
| 10 | 0.797 | 0.520 |
Kesimpulan Simulasi 1: Algoritma mencapai konvergensi pada \(\theta_A \approx 0.796\) dan \(\theta_B \approx 0.520\). Hasil tersebut menunjukkan bahwa Koin A memiliki kecenderungan lebih besar menghasilkan sisi Head, sedangkan Koin B memiliki distribusi yang lebih mendekati koin adil.
Pada simulasi ini digunakan parameter awal dengan perbedaan yang lebih ekstrem, yaitu: \(\theta_A = 0.9\) dan \(\theta_B = 0.1\).
# Inisialisasi Baru
theta_A_new <- 0.9
theta_B_new <- 0.1
log_baru <- data.frame(Iterasi = 0:iterasi, Theta_A = NA, Theta_B = NA)
log_baru[1, 2:3] <- c(theta_A_new, theta_B_new)
for(i in 1:iterasi) {
# E-Step
p_E_ZA <- dbinom(data_H, n_lemparan, theta_A_new)
p_E_ZB <- dbinom(data_H, n_lemparan, theta_B_new)
mu_A <- p_E_ZA / (p_E_ZA + p_E_ZB)
mu_B <- p_E_ZB / (p_E_ZA + p_E_ZB)
# M-Step
theta_A_new <- sum(data_H * mu_A) / sum(mu_A * n_lemparan)
theta_B_new <- sum(data_H * mu_B) / sum(mu_B * n_lemparan)
log_baru[i+1, 2:3] <- c(round(theta_A_new, 3), round(theta_B_new, 3))
}
knitr::kable(log_baru, caption = "Hasil Iterasi dengan Parameter Ekstrem")| Iterasi | Theta_A | Theta_B |
|---|---|---|
| 0 | 0.900 | 0.100 |
| 1 | 0.756 | 0.434 |
| 2 | 0.779 | 0.485 |
| 3 | 0.790 | 0.506 |
| 4 | 0.794 | 0.514 |
| 5 | 0.796 | 0.517 |
| 6 | 0.796 | 0.519 |
| 7 | 0.797 | 0.519 |
| 8 | 0.797 | 0.519 |
| 9 | 0.797 | 0.520 |
| 10 | 0.797 | 0.520 |
Kesimpulan Simulasi 2: Meskipun menggunakan parameter awal yang sangat berbeda, algoritma EM tetap menghasilkan nilai parameter akhir yang mendekati hasil pada Simulasi 1. Hal ini menunjukkan bahwa algoritma EM mampu mencapai solusi yang stabil melalui proses iteratif.
Simulasi berikut dilakukan menggunakan nilai awal parameter yang identik, yaitu \(\theta_A = 0.5\) dan \(\theta_B = 0.5\).
theta_A_sym <- 0.5
theta_B_sym <- 0.5
log_sym <- data.frame(Iterasi = 0:iterasi, Theta_A = NA, Theta_B = NA)
log_sym[1, 2:3] <- c(theta_A_sym, theta_B_sym)
for(i in 1:iterasi) {
p_E_ZA <- dbinom(data_H, n_lemparan, theta_A_sym)
p_E_ZB <- dbinom(data_H, n_lemparan, theta_B_sym)
# Pada titik ini mu_A dan mu_B akan selalu 0.5 karena p_E_ZA == p_E_ZB
mu_A <- p_E_ZA / (p_E_ZA + p_E_ZB)
mu_B <- p_E_ZB / (p_E_ZA + p_E_ZB)
theta_A_sym <- sum(data_H * mu_A) / sum(mu_A * n_lemparan)
theta_B_sym <- sum(data_H * mu_B) / sum(mu_B * n_lemparan)
log_sym[i+1, 2:3] <- c(round(theta_A_sym, 3), round(theta_B_sym, 3))
}
knitr::kable(log_sym, caption = "Hasil Iterasi dengan Parameter Identik (Symmetry Trap)")| Iterasi | Theta_A | Theta_B |
|---|---|---|
| 0 | 0.50 | 0.50 |
| 1 | 0.66 | 0.66 |
| 2 | 0.66 | 0.66 |
| 3 | 0.66 | 0.66 |
| 4 | 0.66 | 0.66 |
| 5 | 0.66 | 0.66 |
| 6 | 0.66 | 0.66 |
| 7 | 0.66 | 0.66 |
| 8 | 0.66 | 0.66 |
| 9 | 0.66 | 0.66 |
| 10 | 0.66 | 0.66 |
Kesimpulan Simulasi 3: Pada simulasi ini, algoritma gagal membedakan karakteristik kedua koin karena nilai awal parameter dibuat identik. Akibatnya, kedua parameter terus bergerak menuju nilai rata-rata global data tanpa mampu melakukan pemisahan distribusi secara optimal.
Ketiga skenario simulasi dirancang untuk menunjukkan pengaruh inisialisasi parameter terhadap proses konvergensi algoritma EM.
Konvergensi Algoritma Pada Simulasi 1 dan Simulasi 2, algoritma EM berhasil mencapai nilai parameter akhir yang relatif sama, yaitu \(\theta_A \approx 0.8\) dan \(\theta_B \approx 0.52\), meskipun menggunakan nilai awal yang berbeda. Hasil ini menunjukkan bahwa algoritma EM mampu menemukan solusi yang stabil melalui proses iteratif sehingga cukup robust terhadap variasi inisialisasi parameter awal. Perlu diperhatikan bahwa pertukaran label parameter (label switching) dapat terjadi pada algoritma EM. Kondisi ini tidak memengaruhi hasil estimasi secara statistik karena label A dan B hanya berfungsi sebagai identitas distribusi.
Inferensi Variabel Laten Berdasarkan estimasi parameter akhir, dapat dilakukan inferensi terhadap kemungkinan identitas koin yang digunakan pada setiap set observasi. Set 2, Set 3, dan Set 5 memiliki jumlah kemunculan Head yang tinggi sehingga lebih mungkin berasal dari Koin A. Sebaliknya, Set 1 dan Set 4 cenderung berasal dari Koin B karena proporsi Head lebih rendah. Proses ini menunjukkan kemampuan EM dalam memperkirakan variabel laten yang tidak diamati secara langsung melalui probabilitas pada tahap E-Step.
Fenomena Symmetry Trap Pada Simulasi 3, algoritma mengalami kegagalan dalam membedakan karakteristik kedua koin. Hal ini terjadi karena nilai awal parameter dibuat identik, yaitu \(\theta_A = \theta_B = 0.5\). Akibatnya, probabilitas tanggung jawab (responsibility) yang dihasilkan pada E-Step selalu sama besar untuk setiap observasi. Kondisi tersebut menyebabkan proses M-Step hanya menghitung rata-rata global seluruh data (\(\frac{33 \text{ Head}}{50 \text{ Total}} = 0.66\)) sehingga kedua parameter tetap bernilai sama dan tidak terjadi pemisahan distribusi. Kondisi ini dikenal sebagai symmetry trap atau jebakan simetri.
Implikasi Praktis Hasil simulasi menunjukkan bahwa pemilihan nilai awal parameter memiliki peran penting dalam implementasi algoritma EM. Inisialisasi parameter yang identik sebaiknya dihindari karena dapat menyebabkan algoritma gagal melakukan pemisahan kelompok data. Oleh karena itu, pemberian nilai awal yang berbeda diperlukan agar proses iterasi dapat mengidentifikasi karakteristik masing-masing distribusi secara lebih optimal.