Analisis Data Kategori merupakan salah satu analisis statistika yang banyak digunakan dan dapat diterapkan untuk menganalisis permasalahan biomedik dan sosial yang semakin kompleks dan bersifat multidimensional.Analisis Data Kategori merupakan mata kuliah yang membahas hubungan fungsional antara satu atau lebih variabel kategori dengan satu atau lebih variabel kategori atau variabel kontinu atau gabungan dari keduanya (berupa model-model linier) dalam kasus error yang berdistribusi famili eksponen, atas dasar statistika sebagai bahan acuan bagi para pengguna khususnya dalam bidang biomedik dan sosial.
Memahami Pola dan Hubungan: Mengidentifikasi pola atau tren dalam data kategorik serta hubungan antara variabel kategori.
Pengelompokan Data: Mengelompokkan data ke dalam kategori yang bermakna untuk mempermudah interpretasi dan analisis lebih lanjut.
Pengambilan Keputusan: Menyediakan informasi yang relevan untuk mendukung proses pengambilan keputusan yang lebih tepat dan efektif
Peningkatan Efisiensi dan Produktivitas: Dengan memahami data kategorik, organisasi dapat mengidentifikasi area yang memerlukan perbaikan, sehingga meningkatkan efisiensi operasional.
Identifikasi Peluang Baru: Analisis data kategorik dapat mengungkap peluang baru, seperti segmen pasar potensial atau kebutuhan pelanggan yang belum terpenuhi.
Pengurangan Risiko: Dengan memahami hubungan antar variabel kategori, organisasi dapat mengidentifikasi potensi risiko dan mengambil langkah pencegahan yang diperlukan
Dalam analisis data kategori, desain sampling memegang peran penting dalam memastikan validitas dan reliabilitas hasil penelitian. Pemilihan desain sampling yang sesuai bergantung pada tujuan penelitian serta jenis data yang dikumpulkan. Secara umum, terdapat dua pendekatan utama dalam desain sampling untuk analisis data kategori, yaitu prospective sampling dan retrospective sampling. Setiap pendekatan memiliki karakteristik serta metode pengambilan sampel yang berbeda dan sering diterapkan dalam penelitian eksperimental maupun observasional, seperti eksperimen, studi kohort, dan studi kasus-kontrol.
Prospective sampling adalah metode pengambilan sampel di mana peneliti mengamati individu atau kelompok tertentu dalam jangka waktu ke depan. Sampel diambil sebelum kejadian atau hasil yang diamati terjadi.
Eksperimen dalam prospective sampling melibatkan intervensi atau perlakuan tertentu terhadap subjek penelitian untuk mengamati dampaknya. Peneliti biasanya membagi sampel ke dalam kelompok eksperimen dan kontrol untuk membandingkan hasilnya.
Studi kohort adalah metode penelitian observasional di mana sekelompok individu yang memiliki karakteristik tertentu diikuti dalam jangka waktu tertentu untuk mengamati perkembangan suatu kondisi atau penyakit. Kohort dapat dibagi menjadi kelompok yang terpapar dan tidak terpapar suatu faktor risiko untuk melihat perbedaannya.
Retrospective sampling adalah metode pengambilan sampel yang dilakukan setelah suatu kejadian telah terjadi. Data yang digunakan biasanya berasal dari catatan sejarah, rekam medis, atau wawancara retrospektif.
Studi kasus-kontrol adalah penelitian observasional di mana dua kelompok dibandingkan: satu kelompok yang mengalami suatu kejadian (kasus) dan satu kelompok yang tidak (kontrol). Peneliti kemudian mencari faktor risiko atau penyebab yang mungkin berkontribusi terhadap perbedaan tersebut.
Studi kohort retrospektif mirip dengan studi kohort prospektif, tetapi dilakukan dengan menggunakan data masa lalu. Peneliti menelusuri kembali data dari sekelompok individu untuk melihat hubungan antara faktor risiko dan hasil yang diamati.
Tabel kontingensi 2x2 adalah tabel yang digunakan untuk menunjukkan hubungan antara dua variabel kategori, masing-masing dengan dua kategori.Tabel ini berbentuk matriks 2 baris dan 2 kolom, di mana sel-selnya berisi jumlah data untuk setiap kombinasi kategori dari kedua variabel tersebut.
Contoh Kasus Misalkan kita memiliki data mengenai hubungan antara aktivitas fisik dan kejadian obesitas.
\[ \begin{array}{|l|c|c|c|} \hline & \textbf{Obesitas (+)} & \textbf{Obesitas (-)} & \textbf{Total} \\ \hline \textbf{Tidak Aktif} & 60 & 40 & 100 \\ \textbf{Aktif} & 30 & 70 & 100 \\ \hline \textbf{Total} & 90 & 110 & 200 \\ \hline \end{array} \]
# Data dalam tabel kontingensi
activity <- matrix(c(60, 40, 30, 70), nrow = 2, byrow = TRUE,
dimnames = list(PhysicalActivity = c("Inactive", "Active"),
Obesity = c("Yes", "No")))
activity
## Obesity
## PhysicalActivity Yes No
## Inactive 60 40
## Active 30 70
\[ P(A \cap B) = \frac{n(A \cap B)}{n(S)} \]
# Peluang bersama
joint_prob <- activity / sum(activity)
joint_prob
## Obesity
## PhysicalActivity Yes No
## Inactive 0.30 0.20
## Active 0.15 0.35
Interpretasi: Peluang bersama menunjukkan proporsi individu dalam setiap kombinasi kategori (aktivitas fisik dan obesitas) terhadap total populasi.
\[ P(A) = \frac{n(A)}{n(S)} \]
# Peluang marginal
marginal_activity <- rowSums(activity) / sum(activity)
marginal_obesity <- colSums(activity) / sum(activity)
list(marginal_activity = marginal_activity, marginal_obesity = marginal_obesity)
## $marginal_activity
## Inactive Active
## 0.5 0.5
##
## $marginal_obesity
## Yes No
## 0.45 0.55
Interpretasi: Peluang marginal menggambarkan proporsi individu yang aktif/tidak aktif dan proporsi individu yang obesitas/tidak obesitas dalam populasi.
\[ P(A | B) = \frac{P(A \cap B)}{P(B)} \]
# Peluang bersyarat berdasarkan aktivitas fisik
conditional_prob_activity <- activity / rowSums(activity)
conditional_prob_activity
## Obesity
## PhysicalActivity Yes No
## Inactive 0.6 0.4
## Active 0.3 0.7
Interpretasi: Peluang bersyarat menunjukkan kemungkinan seseorang obesitas berdasarkan tingkat aktivitas fisiknya.
\[ RD = P(O | I^-) - P(O | I^+) \]
# Beda peluang
RD <- conditional_prob_activity[1,1] - conditional_prob_activity[2,1]
RD
## [1] 0.3
\[ RD = 0.60 - 0.30 = 0.30 \] Interpretasi: Perbedaan peluang menunjukkan bahwa individu yang tidak aktif memiliki 30% lebih besar kemungkinan mengalami obesitas dibandingkan individu yang aktif.
\[ RR = \frac{P(O | I^-)}{P(O | I^+)} \]
# Risiko relatif
RR <- conditional_prob_activity[1,1] / conditional_prob_activity[2,1]
RR
## [1] 2
\[ RR = \frac{0.60}{0.30} = 2.0 \] Interpretasi: Risiko relatif menunjukkan bahwa individu yang tidak aktif memiliki 2 kali lipat risiko obesitas dibandingkan individu yang aktif.
\[ OR = \frac{P(O | I^-) / P(\neg O | I^-)}{P(O | I^+) / P(\neg O | I^+)} \]
# Rasio odds
OR <- (conditional_prob_activity[1,1] / conditional_prob_activity[1,2]) /
(conditional_prob_activity[2,1] / conditional_prob_activity[2,2])
OR
## [1] 3.5
\[ OR = \frac{0.60 / 0.40}{0.30 / 0.70} = \frac{1.5}{0.43} = 3.49 \] Interpretasi: Rasio odds menunjukkan bahwa individu yang tidak aktif memiliki 3.49 kali lipat peluang lebih besar mengalami obesitas dibandingkan individu yang aktif.
Kesimpulan Berdasarkan perhitungan yang telah dilakukan: - Beda peluang (RD) menunjukkan bahwa kurangnya aktivitas fisik meningkatkan kemungkinan obesitas sebesar 30%. - Risiko relatif (RR) menunjukkan bahwa individu yang tidak aktif memiliki risiko 2 kali lipat lebih tinggi terkena obesitas dibandingkan individu yang aktif. - Rasio odds (OR) menunjukkan bahwa individu yang kurang aktif memiliki peluang 3.49 kali lipat lebih besar untuk mengalami obesitas dibandingkan individu yang aktif.
Hasil ini menunjukkan pentingnya aktivitas fisik dalam mengurangi risiko obesitas.
Perbandingan RD, RR, OR
Ukuran | Definisi | Desain Sampling yang Cocok | Interpretasi |
---|---|---|---|
RD (Risk Difference) | Selisih peluang kejadian antara dua kelompok | Studi kohort, studi eksperimental | Mengukur peningkatan atau penurunan risiko absolut |
RR (Relative Risk) | Perbandingan peluang kejadian antara dua kelompok | Studi kohort, studi eksperimental | Mengukur seberapa besar risiko relatif dalam satu kelompok dibandingkan kelompok lain |
OR (Odds Ratio) | Perbandingan odds kejadian antara dua kelompok | Studi kasus-kontrol, studi kohort | Mengukur kekuatan asosiasi antara variabel paparan dan kejadian |
Inferensi dalam statistik merujuk pada proses menarik kesimpulan tentang suatu populasi berdasarkan data sampel. Dalam konteks tabel kontingensi dua arah, inferensi digunakan untuk mengevaluasi hubungan antara dua variabel kategorikal yang disajikan dalam bentuk tabel kontingensi.
Tabel kontingensi sendiri merupakan tabel yang menunjukkan distribusi frekuensi dari dua variabel kategorikal dalam format matriks. Tujuan utama dari inferensi pada tabel ini adalah untuk mengkaji keterkaitan antara kedua variabel melalui estimasi parameter serta pengujian hipotesis.
Inferensi tabel kontingensi dua arah digunakan untuk menilai hubungan antara dua variabel kategori dengan melakukan estimasi dan pengujian hipotesis.
Inferensi dalam tabel kontingensi dua arah dapat dibagi menjadi dua kategori utama:
Estimasi
Pengujian
Estimasi bertujuan untuk memperkirakan parameter populasi berdasarkan data sampel. Estimasi dibagi menjadi :
Estimasi titik adalah estimasi satu nilai parameter populasi berdasarkan sampel yang diperoleh.
\[ \hat{p} = \frac{x}{n} \]
Dimana: - \(\hat{p}\) adalah estimasi titik dari proporsi kejadian - \(x\) adalah jumlah kejadian dalam sampel - \(n\) adalah ukuran sampel
Estimasi interval memberikan rentang kemungkinan nilai parameter populasi berdasarkan data sampel dengan tingkat kepercayaan tertentu.
\[ CI = \hat{p} \pm Z_{\alpha/2} \sqrt{\frac{\hat{p} (1 - \hat{p})}{n}} \]
Dimana: - \(Z_{\alpha/2}\) adalah nilai kritis dari distribusi normal - \(CI\) adalah interval kepercayaan
Uji proporsi digunakan untuk membandingkan proporsi suatu kejadian di antara dua kelompok dalam tabel kontingensi, dengan tujuan mengetahui apakah terdapat perbedaan proporsi kejadian yang signifikan antara kedua kelompok tersebut.
Tabel kontingensi \(2 \times 2\) memiliki struktur sebagai berikut:
\[ \begin{array}{|c|c|c|c|} \hline & \text{Kejadian (+)} & \text{Tidak Kejadian (-)} & \text{Total} \\ \hline \text{Grup 1} & n_{11} & n_{12} & n_{1.} \\ \text{Grup 2} & n_{21} & n_{22} & n_{2.} \\ \hline \text{Total} & n_{.1} & n_{.2} & n \\ \hline \end{array} \]
Untuk menguji hipotesis bahwa tidak ada perbedaan proporsi antara dua kelompok, kita menggunakan uji z dua proporsi, dengan hipotesis:
Estimasi proporsi dalam masing-masing kelompok diberikan oleh:
\[ \hat{p}_1 = \frac{n_{11}}{n_{1.}}, \quad \hat{p}_2 = \frac{n_{21}}{n_{2.}} \]
Estimasi proporsi gabungan (pooling proportion):
\[ \hat{p} = \frac{n_{11} + n_{21}}{n_{1.} + n_{2.}} \]
Statistik uji untuk uji proporsi dua sampel:
\[ Z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{\hat{p}(1 - \hat{p}) \left( \frac{1}{n_{1.}} + \frac{1}{n_{2.}} \right)}} \]
Statistik uji \(Z\) mengikuti distribusi normal baku \(N(0,1)\), dan p-value dihitung berdasarkan nilai kritis dari distribusi normal.
Jika \(|Z|\) lebih besar dari nilai kritis untuk tingkat signifikansi tertentu \(\alpha\) (misalnya 1.96 untuk \(\alpha = 0.05\)), maka hipotesis nol ditolak, yang berarti ada perbedaan signifikan antara dua proporsi.
Uji ini cocok digunakan dalam studi kohort dan eksperimen klinis.
Contoh Kasus
Sebuah penelitian klinis ingin menguji efektivitas obat baru terhadap penyakit tertentu. Dua kelompok pasien diuji:
Hasilnya sebagai berikut:
\[ \begin{array}{|l|c|c|c|} \hline \textbf{Kelompok} & \textbf{Sembuh} & \textbf{Tidak Sembuh} & \textbf{Total} \\ \hline \text{A} & 42 & 18 & 60 \\ \text{B} & 30 & 30 & 60 \\ \hline \textbf{Total} & 72 & 48 & 120 \\ \hline \end{array} \]
Langkah-langkah Perhitungan
Proporsi sampel masing-masing kelompok:
\[ \hat{p}_1 = \frac{\text{Jumlah sembuh pada Obat Baru}}{\text{Total Obat Baru}} = \frac{42}{60} = 0.7 \]
\[ \hat{p}_2 = \frac{\text{Jumlah sembuh pada Plasebo}}{\text{Total Plasebo}} = \frac{30}{60} = 0.5 \]
Jadi,
- Proporsi sembuh pada kelompok Obat Baru: 0.7
- Proporsi sembuh pada kelompok Plasebo: 0.5
Proporsi gabungan seluruh sampel:
\[ \hat{p} = \frac{\text{Jumlah total sembuh}}{\text{Total seluruh pasien}} = \frac{72}{120} = 0.6 \]
Proporsi gabungan sembuh adalah 0.6.
Rumus statistik uji z untuk dua proporsi:
\[ z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{\hat{p}(1 - \hat{p}) \left( \frac{1}{n_1} + \frac{1}{n_2} \right)}} \]
Dengan:
Langkah perhitungan:
\[ \hat{p}(1 - \hat{p}) = 0.6 \times 0.4 = 0.24 \]
\[ \frac{1}{n_1} + \frac{1}{n_2} = \frac{1}{60} + \frac{1}{60} = \frac{2}{60} = \frac{1}{30} \]
\[ \sqrt{0.24 \times \frac{1}{30}} = \sqrt{0.008} = 0.08944 \]
Lalu, hitung nilai z:
\[ z = \frac{0.7 - 0.5}{0.08944} = \frac{0.2}{0.08944} \approx 2.237 \]
Misalkan kita gunakan tingkat signifikansi \(\alpha = 0.05\) untuk uji dua sisi.
Karena:
\[ |z| = 2.237 > 1.96 \]
Maka keputusan: Tolak H0.
Berdasarkan hasil uji proporsi:
Kesimpulan: Ada cukup bukti statistik untuk menyimpulkan bahwa obat baru lebih efektif dibandingkan plasebo dalam meningkatkan kesembuhan pasien.
# Data jumlah sembuh dan total pasien masing-masing kelompok
sembuh <- c(42, 30)
total <- c(60, 60)
# Buat tabel kontingensi
tabel <- matrix(c(42, 18, 30, 30),
nrow = 2,
byrow = TRUE,
dimnames = list(Kelompok = c("Obat Baru", "Plasebo"),
Hasil = c("Sembuh", "Tidak Sembuh")))
# Tampilkan tabel kontingensi
tabel
## Hasil
## Kelompok Sembuh Tidak Sembuh
## Obat Baru 42 18
## Plasebo 30 30
# Hitung proporsi masing-masing kelompok
proporsi <- sembuh / total
proporsi
## [1] 0.7 0.5
# Hitung proporsi gabungan
prop_gabungan <- sum(sembuh) / sum(total)
prop_gabungan
## [1] 0.6
# Manual hitung standar error
SE <- sqrt(prop_gabungan * (1 - prop_gabungan) * (1/total[1] + 1/total[2]))
SE
## [1] 0.08944272
# Hitung statistik z manual di R
z_manual <- (proporsi[1] - proporsi[2]) / SE
z_manual
## [1] 2.236068
# Lakukan uji proporsi dengan prop.test
uji_proporsi <- prop.test(x = sembuh, n = total, correct = FALSE)
uji_proporsi
##
## 2-sample test for equality of proportions without continuity correction
##
## data: sembuh out of total
## X-squared = 5, df = 1, p-value = 0.02535
## alternative hypothesis: two.sided
## 95 percent confidence interval:
## 0.02838653 0.37161347
## sample estimates:
## prop 1 prop 2
## 0.7 0.5
# Ambil nilai p-value dari output
uji_proporsi$p.value
## [1] 0.02534732
Kelompok Obat Baru memiliki tingkat kesembuhan secara signifikan lebih tinggi dibandingkan kelompok Plasebo.
Nilai p-value yang kecil (< 0.05) mengindikasikan bahwa perbedaan proporsi kesembuhan tidak disebabkan oleh kebetulan semata.
Dengan demikian, penelitian ini memberikan bukti statistik bahwa obat baru efektif meningkatkan kesembuhan pasien
Uji asosiasi pada tabel kontingensi 2x2 digunakan untuk mengevaluasi hubungan antara dua variabel kategori. Terdapat tiga ukuran utama yang digunakan dalam uji ini yaitu:
Hipotesis Uji dalam Tabel Kontingensi 2x2
Untuk setiap uji asosiasi, hipotesus yang diuji adalah:
Risk Difference mengukur perbedaan absolut dalam probabilitas kejadian antara dua kelompok:
\[ RD = \left( \frac{n_{11}}{n_{1.}} \right) - \left( \frac{n_{21}}{n_{2.}} \right) \]
Standard Error untuk RD:
\[ SE(RD) = \sqrt{ \frac{ \hat{p}_1 (1 - \hat{p}_1) }{n_{1.}} + \frac{ \hat{p}_2 (1 - \hat{p}_2) }{n_{2.}} } \]
Statistik uji Z:
\[ Z_{RD} = \frac{RD}{SE(RD)} \]
Relative Risk membandingkan kemungkinan kejadian antara dua kelompok:
\[ RR = \frac{ \frac{n_{11}}{n_{1.}} }{ \frac{n_{21}}{n_{2.}} } \]
Standard Error untuk log(RR):
\[ SE(\ln RR) = \sqrt{ \frac{1}{n_{11}} - \frac{1}{n_{1.}} + \frac{1}{n_{21}} - \frac{1}{n_{2.}} } \]
Statistik uji Z:
\[ Z_{RR} = \frac{ \ln RR }{ SE(\ln RR) } \]
Odds Ratio membandingkan peluang kejadian antara dua kelompok:
\[ OR = \frac{n_{11} \times n_{22}}{n_{12} \times n_{21}} \]
Standard Error untuk log(OR):
\[ SE(\ln OR) = \sqrt{ \frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}} } \]
Statistik uji Z:
\[ Z_{OR} = \frac{ \ln OR }{ SE(\ln OR) } \]
Contoh Kasus Misalkan kita ingin mengetahui hubungan antara kebiasaan berolahraga dan kejadian obesitas. Penelitian dilakukan terhadap 200 orang yang dibagi dalam dua kelompok berdasarkan kebiasaan olahraga:
Grup 1: Orang yang suka berolahraga
Grup 2: Orang yang tidak suka berolahraga
\[ \begin{array}{|l|c|c|c|} \hline \textbf{Grup} & \textbf{Obesitas (+)} & \textbf{Obesitas (-)} & \textbf{Total} \\ \hline \text{Grup 1} & 30 & 70 & 100 \\ \text{Grup 2} & 90 & 10 & 100 \\ \hline \textbf{Total} & 120 & 80 & 200 \\ \hline \end{array} \]
Perhitungan Manual untuk RD - Risk Difference
Misalkan \(\hat{p}_1 = \dfrac{30}{100} = 0.3\) dan \(\hat{p}_2 = \dfrac{90}{100} = 0.9\)
Risk Difference: \[RD = 0.3 - 0.9 = -0.6\]
Standard Error : \[ SE(RD) = \sqrt{ \dfrac{0.3(1 - 0.3)}{100} + \dfrac{0.9(1 - 0.9)}{100} } = \sqrt{ \dfrac{0.21}{100} + \dfrac{0.09}{100} } = \sqrt{0.003} = 0.05477 \]
Statistik uji Z: \[Z_{RD} = \dfrac{-0.6}{0.05477} = -10.96\]
Perhitungan Manual untuk RR - Relative Risk: \[RR = \dfrac{0.3}{0.9} = 0.333\]
Standard Error untuk ln RR: \[ SE(\ln RR) = \sqrt{ \dfrac{1}{30} - \dfrac{1}{100} + \dfrac{1}{90} - \dfrac{1}{100} } = \sqrt{0.0333 + 0.0111 - 0.02} = \sqrt{0.0244} = 0.1562\]
Statistik uji Z: \[Z_{RR} = \dfrac{\ln(0.333)}{0.1562} = \dfrac{-1.0986}{0.1562} = -7.033\]
Perhitungan Manual untuk OR - Odds Ratio: \[OR = \dfrac{30 \times 10}{70 \times 90} = \dfrac{300}{6300} = 0.0476\]
Standard Error untuk ln OR \[ SE(\ln OR) = \sqrt{ \dfrac{1}{30} + \dfrac{1}{70} + \dfrac{1}{90} + \dfrac{1}{10} } = \sqrt{0.0333 + 0.0143 + 0.0111 + 0.1} = \sqrt{0.1587} = 0.3983\]
Statistik uji Z
\[ Z_{OR} = \dfrac{\ln(0.0476)}{0.3983} = \dfrac{-3.048}{0.3983} = -7.65\]
Perhitungan menggunakan software R
# Data:
# Grup 1 (suka olahraga): 30 obesitas, 70 tidak obesitas
# Grup 2 (tidak suka olahraga): 90 obesitas, 10 tidak obesitas
n11 <- 30 # suka olahraga dan obesitas
n12 <- 70 # suka olahraga dan tidak obesitas
n21 <- 90 # tidak suka olahraga dan obesitas
n22 <- 10 # tidak suka olahraga dan tidak obesitas
n1. <- n11 + n12 # total grup 1
n2. <- n21 + n22 # total grup 2
# Risk Difference
p1 <- n11 / n1.
p2 <- n21 / n2.
rd <- p1 - p2
se_rd <- sqrt((p1 * (1 - p1) / n1.) + (p2 * (1 - p2) / n2.))
z_rd <- rd / se_rd
# Relative Risk
rr <- p1 / p2
se_ln_rr <- sqrt((1/n11) - (1/n1.) + (1/n21) - (1/n2.))
z_rr <- log(rr) / se_ln_rr
# Odds Ratio
or <- (n11 * n22) / (n12 * n21)
se_ln_or <- sqrt(1/n11 + 1/n12 + 1/n21 + 1/n22)
z_or <- log(or) / se_ln_or
# Output
list(
RD = rd,
SE_RD = se_rd,
Z_RD = z_rd,
RR = rr,
SE_Ln_RR = se_ln_rr,
Z_RR = z_rr,
OR = or,
SE_Ln_OR = se_ln_or,
Z_OR = z_or
)
## $RD
## [1] -0.6
##
## $SE_RD
## [1] 0.05477226
##
## $Z_RD
## [1] -10.95445
##
## $RR
## [1] 0.3333333
##
## $SE_Ln_RR
## [1] 0.1563472
##
## $Z_RR
## [1] -7.026748
##
## $OR
## [1] 0.04761905
##
## $SE_Ln_OR
## [1] 0.3984095
##
## $Z_OR
## [1] -7.641691
Kesimpulan
Risk Difference (RD) mengukur perbedaan risiko absolut. Dalam kasus ini, RD = -0,60, yang berarti risiko obesitas pada kelompok yang suka olahraga lebih rendah sebesar 60% dibandingkan kelompok yang tidak suka olahraga.
Relative Risk (RR) membandingkan kemungkinan kejadian antara dua kelompok. Hasil perhitungan menunjukkan RR = 0,33, artinya individu yang suka olahraga memiliki sepertiga risiko obesitas dibandingkan mereka yang tidak suka olahraga.
Odds Ratio (OR) membandingkan peluang kejadian antara dua kelompok. Diperoleh OR = 0,0476, yang berarti peluang obesitas pada orang yang suka olahraga hanya sekitar 4,76% dari peluang obesitas pada orang yang tidak suka olahraga.
Standard error dan statistik uji Z digunakan untuk menilai signifikansi statistik masing-masing ukuran asosiasi. Semua nilai Z yang diperoleh untuk RD, RR, dan OR lebih besar dari 6 dalam nilai mutlak (Z_RD = -7,75; Z_RR = -6,85; Z_OR = -6,97), yang berarti seluruh ukuran asosiasi tersebut signifikan secara statistik (p < 0,001).
Uji independensi digunakan untuk menentukan apakah ada hubungan statistik antara dua variabel kategorikal.
Uji Chi-Square digunakan untuk menguji apakah ada hubungan antara dua variabel kategorikal
Rumus Chi-Square : \[X^2 = \sum \frac{(O - E)^2}{E}\]
di mana : - \(O\) adalah nilai
observasi dalam tabel kontingensi.
- \(E\) adalah nilai yang diharapkan,
dihitung sebagai:
\[E_{ij} = \frac{R_i \times C_j}{N}\]
dengan
Contoh Kasus
Dalam sebuah studi mengenai kesehatan mental mahasiswa, peneliti ingin mengetahui apakah ada hubungan antara tingkat stres dan kesehatan mental. Penelitian ini melibatkan 100 mahasiswa dari berbagai jurusan di sebuah universitas.
Peneliti membagi mahasiswa ke dalam dua kategori berdasarkan tingkat stres(Stres Tinggi dan Stres Rendah) Kemudian, peneliti juga mencatat kondisi kesehatan mental mahasiswa, diklasifikasikan menjadi (Sehat dan Mengalami Gangguan Mental
Peneliti ingin mengetahui apakah tingkat stres berasosiasi dengan gangguan kesehatan mental.
\[ \begin{array}{|l|c|c|c|} \hline \textbf{Tingkat Stres} & \textbf{Sehat} & \textbf{Gangguan Mental} & \textbf{Total} \\ \hline \text{Stres Tinggi} & 20 & 40 & 60 \\ \text{Stres Rendah} & 30 & 10 & 40 \\ \hline \textbf{Total} & 50 & 50 & 100 \\ \hline \end{array} \]
Hitung nilai yang diharapkan \(E\): \[E_{11} = \frac{60 \times 50}{100} = 30 \quad\] \[E_{12} = \frac{60 \times 50}{100} = 30\] \[E_{21} = \frac{40 \times 50}{100} = 20 \quad\] \[E_{22} = \frac{40 \times 50}{100} = 20 \]
Lalu, hitung \(X^2\): \[ \chi^2 = \frac{(20 - 30)^2}{30} + \frac{(40 - 30)^2}{30} + \frac{(30 - 20)^2}{20} + \frac{(10 - 20)^2}{20} \]
\[ = \frac{100}{30} + \frac{100}{30} + \frac{100}{20} + \frac{100}{20} \]
\[ = 3.33 + 3.33 + 5 + 5 = 16.66 \] Dengan derajat kebebasan \(df = (2-1)(2-1) = 1\), kita membandingkan dengan tabel \(\chi^2\), yaitu \(3.841\) pada \(\alpha = 0.05\).
Karena \(16.66 > 3.841\), maka tolak hipotesis nol. Artinya Terdapat hubungan yang signifikan antara tingkat stres dan kesehatan mental mahasiswa.
set.seed(123)
data <- matrix(c(20, 40, 30, 10), nrow = 2, byrow = TRUE)
dimnames(data) <- list( "Tingkat_Stres" = c("Tinggi", "Rendah"),"Kesehatan_Mental" = c("Sehat", "Gangguan"))
print(data)
## Kesehatan_Mental
## Tingkat_Stres Sehat Gangguan
## Tinggi 20 40
## Rendah 30 10
# Uji Chi-Square
chisq_test <- chisq.test(data)
print(chisq_test)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: data
## X-squared = 15.042, df = 1, p-value = 0.0001052
Dengan perhitungan menggunakan software R diperoleh nilai pvalue sebesar 0,0001052 lebih kecil dari alpha (0,05). Artinya, terdapat hubungan antara variabel tingkat stress dan kesehatan mahasiswa.
Partisi Chi-Square digunakan untuk mengetahui kategori mana dalam tabel kontingensi yang memiliki andil terhadap hubungan yang signifikan. Ketika uji Chi-Square pada tabel kontingensi berukuran I × J menunjukkan hasil signifikan, maka partisi Chi-Square memungkinkan kita untuk menguraikan pengaruh tersebut ke dalam subkelompok yang lebih kecil.
Menurut Lancaster dan Irwin, statistik uji chi-square dapat diuraikan menjadi beberapa komponen sesuai dengan derajat bebas dalam tabel kontingensi. Ini memungkinkan kita mengidentifikasi kategori-kategori yang berkontribusi secara signifikan terhadap hubungan tersebut.
Fenomena ini juga berkaitan dengan Paradoks Simpson, yaitu situasi ketika tren yang tampak dalam beberapa kelompok data bisa menghilang atau bahkan berbalik arah saat data tersebut digabungkan.
Langkah-langkah dalam Partisi Chi-Square:
Gunakan tabel I × J dan hitung nilai Chi-Square total.
Pecah tabel tersebut menjadi beberapa tabel 2 × 2, sebanyak (I-1)(J-1) buah.
Hitung nilai Chi-Square untuk masing-masing tabel 2 × 2.
Interpretasikan kategori mana yang memberikan kontribusi paling signifikan.
Contoh kasus Sebuah lembaga survei pendidikan ingin mengetahui apakah ada hubungan antara Jenis Kelamin dan Tipe Aktivitas Ekstrakurikuler yang Dipilih oleh mahasiswa diuniversitas-universitas Indonesia. Mereka mensurvei 5000 mahasiswa dari berbagai jurusan, dan menanyakan jenis kegiatan ekstrakurikuler mana yang paling mereka ikuti atau minati: Organisasi Kemahasiswaan, Kegiatan Sosial/Relawan, atau Tidak Aktif dalam Kegiatan Ekstrakurikuler.
Peneliti ingin mengetahui apakah keterlibatan dalam aktivitas di luar akademik tersebut dipengaruhi oleh jenis kelamin mahasiswa.
\[ \begin{array}{|l|c|c|c|c|} \hline \textbf{Gender} & \textbf{Organisasi} & \textbf{Sosial/Relawan} & \textbf{Tidak aktif} & \textbf{Total} \\ \hline \text{Female} & 920 & 710 & 870 & 2500 \\ \text{Male} & 1030 & 590 & 880 & 2500 \\ \hline \textbf{Total} & 1950 & 1300 & 1750 & 5000 \\ \hline \end{array} \]
Hipotesis yang diuji: - \(H_0\) : Tidak ada hubungan antara jenis kelamin dan tipe aktivitas ekstrakulikuler. - \(H_1\) : Ada hubungan antara jenis kelamin dan tipe aktivitas ekstrakulikuler.
Langkah 1: Melakukan perhitungan chisquare secara keseluruhan Perhitungan manual:
\[E_{ij} = \frac{(n_{i\cdot} \times n_{\cdot j})}{n}\]
Female, Organisasi: \[ E_{11} = \frac{(2500 \times 1950)}{5000} = 975.00 \]
Female, Sosial/Relawan: \[ E_{12} = \frac{(2500 \times 1300)}{5000} = 650.00 \]
Female, Tidak Aktif: \[ E_{13} = \frac{(2500 \times 1750)}{5000} = 875.00 \]
Male, Organisasi: \[ E_{21} = \frac{(2500 \times 1950)}{5000} = 975.00 \]
Male, Sosial/Relawan: \[ E_{22} = \frac{(2500 \times 1300)}{5000} = 650.00 \]
Male, Tidak Aktif: \[ E_{23} = \frac{(2500 \times 1750)}{5000} = 875.00 \]
\[ \chi^2 = \frac{(920 - 975)^2}{975} + \frac{(710 - 650)^2}{650} + \frac{(870 - 875)^2}{875} + \frac{(1030 - 975)^2}{975} + \frac{(590 - 650)^2}{650} + \frac{(880 - 875)^2}{875} \]
\[ \chi^2 = \frac{(-55)^2}{975} + \frac{(60)^2}{650} + \frac{(-5)^2}{875} + \frac{(55)^2}{975} + \frac{(-60)^2}{650} + \frac{(5)^2}{875} \]
\[ \chi^2 = \frac{3025}{975} + \frac{3600}{650} + \frac{25}{875} + \frac{3025}{975} + \frac{3600}{650} + \frac{25}{875} \]
\[ \chi^2 = 3.10 + 5.54 + 0.03 + 3.10 + 5.54 + 0.03 = \boxed{17.34} \]`
\[\text{Derajat bebas (df)} = (2-1)(3-1) = 2 \] \[\text{Nilai kritis } \chi^2_{\alpha=0.05, df=2} = 5.99\]
Karena \(\chi^2 = 17.34 > 5.99\), maka kita tolak hipotesis nol.Oleh karena itu, terdapat hubunga antara Jenis Kelamin dan Tipe Aktivitas Ekstrakulikuler.
Perhitungan menggunakan software R
#data
data <- matrix(c(920, 710, 870, 1030, 590, 880), nrow = 2, byrow = TRUE)
colnames(data) <- c("Organisasi", "Sosial/Relawan","Tidak Aktif")
rownames(data) <- c("Female", "Male")
print(data)
## Organisasi Sosial/Relawan Tidak Aktif
## Female 920 710 870
## Male 1030 590 880
# Uji Chi-Square
chi_test <- chisq.test(data)
chi_test
##
## Pearson's Chi-squared test
##
## data: data
## X-squared = 17.339, df = 2, p-value = 0.0001717
Dengan perhitungan menggunakan software R diperoleh nilai pvalue sebesar 0,0001717 lebih kecil dari alpha (0,05). Artinya, terdapat hubungan antara Jenis Kelamin dan Tipe Aktivitas Ekstrakulikuler.
Langkah 2: Lakukan partisi chi-square
Partisi 1 : Organisasi vs Sosial/Relawan
\[ \begin{array}{|l|c|c|c|} \hline \textbf{Gender} & \textbf{Organisasi} & \textbf{Sosial/Relawan} & \textbf{Total} \\ \hline \text{Female} & 920 & 710 & 1630 \\ \text{Male} & 1030 & 590 & 1620 \\ \hline \textbf{Total} & 1950 & 1300 & 3250 \\ \hline \end{array} \]
Chi-Square : \(X^2 = 17,52\)
Partisi 2 : Aktif ( Organisasi dan Sosial/Relawan) vs Tidak Aktif\[ \begin{array}{|l|c|c|c|} \hline \textbf{Gender} & \textbf{ Aktif } & \textbf{Tidak Aktif} & \textbf{Total} \\ \hline \text{Female} & 1630 & 870 & 2500 \\ \text{Male} & 1620 & 880 & 2500 \\ \hline \textbf{Total} & 3250 & 1750 & 500 \\ \hline \end{array} \]
Chi-Square : \(X^2 = 0,088\)
Perhitungan menggunakan software R
#data parrisi 1
data_partisi1 <- matrix(c(920, 710, 1030, 590), nrow = 2, byrow = TRUE)
colnames(data_partisi1) <- c("Organisasi", "Sosial/Relawan")
rownames(data_partisi1) <- c("Female", "Male")
print(data_partisi1)
## Organisasi Sosial/Relawan
## Female 920 710
## Male 1030 590
# Uji Chi-Square1
chi_test1 <- chisq.test(data_partisi1)
chi_test1
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: data_partisi1
## X-squared = 16.955, df = 1, p-value = 3.827e-05
#data partisi 2
data_partisi2 <- matrix(c(1630,870,1620,880), nrow = 2, byrow = TRUE)
colnames(data_partisi2) <- c("Aktif", "Tidak Aktif")
rownames(data_partisi2) <- c("Female", "Male")
print(data_partisi2)
## Aktif Tidak Aktif
## Female 1630 870
## Male 1620 880
# Uji Chi-Square2
chi_test2 <- chisq.test(data_partisi2)
chi_test2
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: data_partisi2
## X-squared = 0.071209, df = 1, p-value = 0.7896
Interpretasi Hasil
Hasil uji Chi-Square pertama menunjukkan bahwa terdapat perbedaan yang signifikan antara preferensi partisipasi jenis kegiatan (Organisasi vs Sosial/Relawan) antara mahasiswa laki-laki dan perempuan. Nilai statistik \(X^2 = 16,955\) dengan derajat kebebasan 1, dan nilai p-value sebesar 3.827e-05 (p < 0.001), yang berarti hasil tersebut signifikan secara statistik.Artinya, jenis kelamin memengaruhi pilihan antara aktif dalam organisasi atau kegiatan sosial/relawan.
Hasil uji Chi-Square kedua menunjukkan bahwa tidak terdapat perbedaan yang signifikan dalam proporsi keaktifan (Aktif vs Tidak Aktif) antara mahasiswa laki-laki dan perempuan. Nilai statistik \(X^2 = 0,071209\) dengan derajat kebebasan 1, dan p-value sebesar 0.7896, jauh di atas ambang batas signifikansi umum (α = 0.05).Hal ini menunjukkan bahwa jenis kelamin tidak secara signifikan memengaruhi tingkat keaktifan mahasiswa dalam kegiatan secara keseluruhan.
Kesimpulan Jenis kelamin berperan dalam memengaruhi bentuk kegiatan yang dipilih (antara organisasi atau relawan), namun tidak berpengaruh signifikan terhadap kecenderungan mahasiswa untuk aktif atau tidak aktif dalam kegiatan kemahasiswaan secara umum.
Uji Likelihood Ratio (\(G^2\)) adalah alternatif dari uji chi-square yang digunakan untuk menguji hipotesis independensi dalam tabel kontingensi \(I \times J\). Statistik uji ini diberikan oleh:
\[ G^2 = 2 \sum_i \sum_j n_{ij} \ln \left( \frac{n_{ij}}{\hat{\mu}_{ij}} \right) \] Dimana
\(n_{ij}\) adalah frekuensi
observasi dalam tabel kontingensi.
\(\hat{\mu}_{ij} = n \cdot p_{i+} \cdot p_{+j}\) adalah frekuensi ekspektasi.
Statistik \(G^2\) mengikuti distribusi chi-square dengan derajat bebas(I−1 )(J−1).
Tolak \(H_0\) jika: \[G^2 \geq \chi^2_{(1 - \alpha/2)(I - 1)(J - 1)}\]
Contoh Kasus Sebuah pusat kesehatan masyarakat (puskesmas) ingin meneliti hubungan antara kebiasaan sarapan dan kejadian anemia pada remaja putri usia 15–18 tahun. Data dikumpulkan dari 1000 responden melalui pemeriksaan hemoglobin dan survei pola makan.
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{Anemia (+)} & \textbf{Anemia (-)} & \textbf{Total} \\ \hline \text{Sarapan Teratur} & 120 & 580 & 700 \\ \text{Tidak Sarapan} & 180 & 120 & 300 \\ \hline \textbf{Total} & 300 & 700 & 1000 \\ \hline \end{array} \]
Langkah 1 : Hitung Frekuensi Ekspektasi \[ \hat{\mu}_{11} = \frac{(700 \times 300)}{1000} = 210 \quad \hat{\mu}_{12} = \frac{(700 \times 700)}{1000} = 490 \] \[ \hat{\mu}_{21} = \frac{(300 \times 300)}{1000} = 90 \quad \hat{\mu}_{22} = \frac{(300 \times 700)}{1000} = 210 \]
Langkah 2 : Hitung Statistik Uji \(G^2\)
\[ G^2 = 2 \left[ 120 \ln\left(\frac{120}{210}\right) + 580 \ln\left(\frac{580}{490}\right) + 180 \ln\left(\frac{180}{90}\right) + 120 \ln\left(\frac{120}{210}\right) \right] \]
\[ G^2 = 2 \left[ 120 \ln(0.571) + 580 \ln(1.184) + 180 \ln(2.000) + 120 \ln(0.571) \right] \]
\[ G^2 = 2 \left[ 120 (-0.5606) + 580 (0.169) + 180 (0.6931) + 120 (-0.5606) \right] \]
\[ G^2 = 2 \left[ -67.27 + 98.02 + 124.76 - 67.27 \right] = 2 \times 88.24 = 176.48 \] Langkah 3 : Bandingkan dengan Distribusi Chi-Square Derajat bebas
\[ (I - 1)(J - 1) = (2 - 1)(2 - 1) = 1 \]
Nilai kritis \(\chi^2\) untuk \(\alpha = 0.05\) dengan \(df = 1\) adalah 3.841}.Karena \(G^2 = 176.48 > 3.841\), maka kita tolak hipotesis nol, dan menyimpulkan bahwa ada hubungan signifikan antara kebiasaan sarapan dan kejadian anemia pada remaja putri.
#input data
data <- matrix(c(120, 580, 180, 120), nrow = 2, byrow = TRUE)
colnames(data) <- c("Anemia (+)", "Tidak Anemia (-)")
rownames(data) <- c("Sarapan Teratur", "Tidak Sarapan")
data
## Anemia (+) Tidak Anemia (-)
## Sarapan Teratur 120 580
## Tidak Sarapan 180 120
#hitung frekuensi ekspektasi
data_expected <- chisq.test(data)$expected
#Hitung statistik g^2
G2 <- 2 * sum(data * log(data / data_expected))
#Nilai Kritis chi square untuk df = 1 dan alpha = 0,05
critical_value <- qchisq(0.95, df = 1)
data_expected
## Anemia (+) Tidak Anemia (-)
## Sarapan Teratur 210 490
## Tidak Sarapan 90 210
G2
## [1] 176.5198
critical_value
## [1] 3.841459
Nilai statistik uji \(G^2 =
176.5198\)
Nilai kritis chi-square untuk \(\alpha =
0.05\) dan \(df = 1\) adalah
\(\chi^2_{0.05,1} = 3.841459\)
Karena: \(G^2 = 176.5198 > 3.841459\) Maka, kita menolak H₀
Kesimpulan
Uji Exact Fisher adalah metode statistik yang digunakan untuk menguji independensi antara dua variabel kategorik dalam tabel kontingensi, terutama ketika ukuran sampel kecil. Uji ini disebut “exact” karena tidak bergantung pada pendekatan distribusi asimtotik seperti uji chi-square.
Tujuan utama dari Uji Exact Fisher adalah untuk menentukan apakah terdapat hubungan signifikan antara dua variabel kategorik dalam tabel kontingensi 2×2, terutama ketika:
Kelebihan - Akurat untuk sampel kecil: Tidak seperti chi-square, uji ini memberikan hasil yang valid meskipun ukuran sampel sangat kecil. - Tidak memerlukan asumsi distribusi normal. - Cocok digunakan saat data tidak memenuhi asumsi-asumsi klasik uji statistik parametrik.
Kekurangan - Komputasi intensif: Untuk tabel yang lebih besar dari 2×2, perhitungan menjadi sangat rumit. - Kurang efisien untuk sampel besar: Waktu komputasi meningkat drastis dan menjadi tidak praktis dibandingkan uji lain seperti chi-square.
Distribusi Hipergeometrik Distribusi Hipergeometrik digunakan untuk menghitung probabilitas mendapatkan sejumlah keberhasilan tertentu dalam n pengambilan dari populasi terbagi dua (sukses vs gagal), tanpa pengembalian.
\[ P(X = x) = \frac{\dbinom{K}{x} \dbinom{N-K}{n-x}}{\dbinom{N}{n}} \]
Contoh kasus
Misalkan kita memiliki sebuah populasi sebanyak 𝑁= 50, yang terdiri atas 𝐾= 34 bola putih dan sisanya adalah bola hitam. Dari populasi tersebut, kita mengambil sebuah sampel sebanyak 𝑛= 20 bola, secara acak tanpa pengembalian. Ternyata, dari 20 bola yang terambil tersebut, sebanyak𝑥= 15 bola merupakan bola putih. Kita tertarik untuk mengetahui peluang terambilnya 15 bola putih dari 20 bola dalam kondisi tersebut, yang mengikuti distribusi hipergeometrik.
Menggunakan rumus di atas:
\[ \begin{aligned} P(X = 15) &= \frac{\dbinom{34}{15} \cdot \dbinom{16}{5}}{\dbinom{50}{20}} \\ &= \frac{9,\!091,\!563,\!160 \times 4,\!368}{230,\!847,\!726,\!611} \\ &= \frac{39,\!723,\!870,\!268,\!480}{230,\!847,\!726,\!611} \\ &\approx \boxed{0.1720} \end{aligned} \]
Dengan menggunakan software R, diperoleh nilai sebagai berikut:
#definisi parameter
N <- 50
K <- 34
n <- 20
x <- 15
#hitung probabilitas P(X=15)
dhyper(x, m = K, n = N - K, k = n)
## [1] 0.1720136
Misalkan kita memiliki data berikut: \[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 15 & 5 & 20 \\ \text{Group B} & 19 & 11 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
Langkah perhitungan manual : 1. Asumsikan hipotesis non (tidak ada hubungan antara dua varaibel.) 2. Hitung probabilitas tabel yang diamati menggunakan distribusi hipergeometrik:
choose(34,15)*choose(16,5)/choose(50,20)
## [1] 0.1720136
Berikut semua kombinasi tabel 2x2 yang mungkin dengan total baris dan kolom sama :
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 14 & 6 & 20 \\ \text{Group B} & 20 & 10 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,14)*choose(16,6)/choose(50,20)
## [1] 0.2365187
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 13 & 7 & 20 \\ \text{Group B} & 21 & 9 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,13)*choose(16,7)/choose(50,20)
## [1] 0.2252559
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 12 & 8 & 20 \\ \text{Group B} & 22 & 8 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,12)*choose(16,8)/choose(50,20)
## [1] 0.149744
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 11 & 9 & 20 \\ \text{Group B} & 23 & 7 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,11)*choose(16,9)/choose(50,20)
## [1] 0.06944649
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 10 & 10 & 20 \\ \text{Group B} & 24 & 6 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,10)*choose(16,10)/choose(50,20)
## [1] 0.02228075
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 9 & 11 & 20 \\ \text{Group B} & 25 & 5 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,9)*choose(16,11)/choose(50,20)
## [1] 0.004861254
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 8 & 12 & 20 \\ \text{Group B} & 26 & 4 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,8)*choose(16,12)/choose(50,20)
## [1] 0.0007011424
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 7 & 13 & 20 \\ \text{Group B} & 27 & 3 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,7)*choose(16,13)/choose(50,20)
## [1] 6.392182e-05
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 6 & 14 & 20 \\ \text{Group B} & 28 & 2 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,6)*choose(16,14)/choose(50,20)
## [1] 3.424383e-06
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 5 & 15 & 20 \\ \text{Group B} & 29 & 1 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,5)*choose(16,15)/choose(50,20)
## [1] 9.446574e-08
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{ Ya } & \textbf{Tidak} & \textbf{Total} \\ \hline \text{Group A} & 4 & 16 & 20 \\ \text{Group B} & 30 & 0 & 30 \\ \hline \textbf{Total} & 34 & 16 & 50 \\ \hline \end{array} \]
choose(34,4)*choose(16,16)/choose(50,20)
## [1] 9.840181e-10
Menggunakan perhitungan fisher test dengan software R
data <- matrix(c(15,5,19,11), nrow =2, byrow = TRUE)
data
## [,1] [,2]
## [1,] 15 5
## [2,] 19 11
fisher.test(data)
##
## Fisher's Exact Test for Count Data
##
## data: data
## p-value = 0.5382
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.4291384 7.7565734
## sample estimates:
## odds ratio
## 1.718029
Berdasarkan hasil uji Fisher yang dilakukan, diperoleh nilai p-value = 0.5382.
Karena nilai p lebih besar dari 0.05 (p > 0.05), maka kita tidak dapat menolak hipotesis nol (\(H_0\)).
Ini berarti bahwa tidak terdapat cukup bukti statistik untuk menyatakan adanya hubungan yang signifikan antara kedua variabel yang diuji.
Dengan demikian, kita menyimpulkan bahwa perbedaan proporsi antara kelompok tidak signifikan secara statistik pada tingkat signifikansi 5%.
Uji | Tujuan | Kriteria Data |
---|---|---|
Uji Chi-Square | Menguji apakah dua variabel kategori saling bebas (independen) | Tabel kontingensi besar (I×J), frekuensi harapan tiap sel ≥ 5 |
Partisi Chi-Square | Mengidentifikasi kontribusi sel/kategori tertentu terhadap hasil uji | Setelah hasil Chi-Square signifikan, digunakan untuk analisis lanjutan |
Uji Likelihood Ratio | Alternatif dari Chi-Square; membandingkan model penuh vs terbatas | Cocok untuk model log-linear dan analisis tabel kontingensi multivariat |
Uji Exact Fisher | Menguji independensi secara eksak dalam tabel kecil (biasanya 2×2) | Frekuensi harapan kecil (<5), distribusi hipergeometrik, tabel 2×2 |
Residual digunakan untuk mengidentifikasi sel mana dalam tabel kontingensi yang paling banyak berkontribusi terhadap adanya hubungan (asosiasi) antara dua variabel kategori.
Secara sederhana, residual adalah selisih antara frekuensi observasi aktual dan frekuensi yang diharapkan (ekspektasi) berdasarkan model independensi (di mana diasumsikan tidak ada hubungan antara variabel).
Residual besar (baik positif maupun negatif) → menunjukkan bahwa frekuensi aktual sangat berbeda dari yang diharapkan. Artinya, ada sel yang sangat menyimpang dari ekspektasi → berpotensi menunjukkan adanya hubungan kuat antara kategori.
Residual kecil atau mendekati nol → menunjukkan bahwa nilai observasi hampir sama dengan nilai ekspektasi → tidak ada kontribusi berarti terhadap hubungan variabel → mendukung asumsi independensi.
Jika Residual = 0 - Tidak ada perbedaan antara data aktual dan ekspektasi. - Menunjukkan tidak ada pola keterkaitan yang berarti antara baris dan kolom tabel → tidak ada hubungan yang signifikan antar variabel kategori
Jika residual positif besar - Frekuensi observasi jauh lebih besar dari yang diharapkan. - Menunjukkan adanya hubungan positif yang kuat antar kategori → kemungkinan adanya asosiasi langsung atau korelasi positif.
Jika residual negatif besar - Frekuensi observasi jauh lebih kecil dari yang diharapkan. - Menunjukkan hubungan negatif atau bahkan tidak adanya hubungan antar kategori tertentu → bisa jadi terjadi penghindaran atau pengelompokan terpisah.
Pearson Residual
Pearson Residual digunakan untuk mengukur seberapa besar penyimpangan nilai observasi (\(O_{ij}\)) dari nilai yang diharapkan (\(E_{ij}\)) berdasarkan model independensi dalam tabel kontingensi. Rumus Pearson Residual sebagai berikut :
\[e_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}}}\]
dengan - \(O_{ij}\) adalah nilai observasi pada sel ke-\(i,j\), - \(E_{ij}\) adalah nilai ekspektasi pada sel ke-\(i,j\) yang dihitung berdasarkan asumsi bahwa variabel bersifat independen.
Residual Pearson digunakan untuk mendeteksi sel mana yang memiliki kontribusi besar terhadap nilai statistik uji \(\chi^2\).
Standardized Residual (Adjusted Residual)
Standardized residual, juga dikenal sebagai Adjusted Residual, adalah bentuk standar dari Pearson residual yang mempertimbangkan efek dari distribusi marginal baris dan kolom. Ini memberikan interpretasi yang lebih akurat karena memperhitungkan variasi tambahan.
\[r_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij} (1 - p_{i+})(1 - p_{+j})}}\]
Residual ini sering digunakan dalam analisis post-hoc setelah uji \(\chi^2\) untuk mengidentifikasi sel mana yang menyimpang secara signifikan. Jika \(|r_{ij}| > 2\), maka dapat dianggap signifikan pada taraf 5%.
Contoh Kasus Sebuah survei kecil dilakukan terhadap 35 mahasiswa untuk mengetahui preferensi mereka dalam memilih metode pembelajaran di masa perkuliahan. Mahasiswa dikelompokkan berdasarkan angkatan (Angkatan 2022 dan Angkatan 2023), dan mereka ditanya apakah lebih menyukai metode Kuliah Luring (Tatap Muka) atau Kuliah Daring (Online).
Hasil wawancara menunjukkan bahwa preferensi metode belajar bisa jadi dipengaruhi oleh pengalaman awal mereka saat masuk kuliah. Mahasiswa angkatan 2022 mengalami perkuliahan awal di masa transisi pandemi, sedangkan angkatan 2023 sudah mulai kuliah di era pasca-pandemi.
\[ \begin{array}{|l|c|c|c|} \hline \textbf{} & \textbf{Kuliah Luring} & \textbf{Kuliah Daring} & \textbf{Total} \\ \hline \text{Angkatan 2022} & 10 & 8 & 18 \\ \text{Angkaran 2023} & 6 & 11 & 17 \\ \hline \textbf{Total} & 16 & 19 & 35 \\ \hline \end{array} \]
Menghitung Pearson Residual \(e_{ij}\)
Menghitung nilai ekspektasi \(E_{ij}\) menggunakan rumus : \[ E_{ij} = \frac{(\text{Total Baris}) \times (\text{Total Kolom})}{\text{Total Keseluruhan}} \]
Hitung setiap sel:
\[ e_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}}} \] 4. Hitung setiap sel:
Semua nilai residual berada di sekitar nol (antara -1 dan 1), yang menunjukkan bahwa tidak ada perbedaan yang signifikan antara observasi dan ekspektasi. Dengan demikian, tidak ada bukti kuat adanya hubungan antara angkatan mahasiswa dan preferensi metode belajar, dan data mendukung asumsi independensi antar kategori.
Perhitungan menggunakan software R
# Data Observasi
data <- matrix(c(10, 8, 6, 11), nrow = 2, byrow = TRUE)
data
## [,1] [,2]
## [1,] 10 8
## [2,] 6 11
# Hitung nilai ekspektasi
expected <- chisq.test(data)$expected
# Pearson Residual
pearson_residual <- (data - expected) / sqrt(expected)
# Gunakan expected untuk menghitung p_row dan p_col
row_sum_exp <- rowSums(expected)
col_sum_exp <- colSums(expected)
total_sum <- sum(expected)
p_row <- row_sum_exp / total_sum
p_col <- col_sum_exp / total_sum
# Buat matriks proporsi marginal
p_row_mat <- matrix(rep(p_row, each = length(p_col)), nrow = 2)
p_col_mat <- matrix(rep(p_col, times = length(p_row)), nrow = 2, byrow = TRUE)
# Standardized Residual
standardized_residual <- (data - expected) / sqrt(expected * (1 - p_row_mat) * (1 - p_col_mat))
# Tampilkan hasil
list(Pearson_Residual = pearson_residual, Standardized_Residual = standardized_residual)
## $Pearson_Residual
## [,1] [,2]
## [1,] 0.6175348 -0.5666888
## [2,] -0.6354380 0.5831179
##
## $Standardized_Residual
## [,1] [,2]
## [1,] 1.202619 -1.168736
## [2,] -1.237485 1.202619
Outlier dalam konteks analisis data kategorik merujuk pada sel-sel dalam tabel kontingensi yang memiliki nilai residual yang sangat besar, baik bernilai positif maupun negatif. Residual ini merupakan selisih antara frekuensi yang diamati (observed frequency) dengan frekuensi yang diharapkan (expected frequency) berdasarkan asumsi independensi antar variabel.
Sel-sel yang dianggap sebagai outlier menunjukkan bahwa terdapat penyimpangan signifikan dari pola yang diprediksi oleh model independensi. Artinya, terdapat kombinasi kategori yang muncul jauh lebih sering atau jauh lebih jarang dibandingkan dengan yang seharusnya terjadi bila tidak ada hubungan antar variabel.
Menggunakan residual untuk mendeteksi outlier
\[e_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}}}\] Jika nilai \(|e_{ij}| > 2\), maka sel tersebut dianggap sebagai indikasi adanya outlier.
\[r_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij} (1 - p_{i+})(1 - p_{+j})}}\]
Jika nilai \(|r_{ij}| > 2\), maka sel tersebut dianggap sebagai outlier signifikan.
Misalkan kita memiliki tabel kontingensi berikut:
\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Respon Positif} & \textbf{Respon Negatif} & \textbf{Total} \\ \hline \textbf{Kelompok 1} & 60 & 30 & 90 \\ \textbf{Kelompok 2} & 20 & 90 & 110 \\ \hline \textbf{Total} & 80 & 120 & 200 \\ \hline \end{array} \]
Setelah menghitung residual, kita menemukan:
\[ \begin{array}{|c|c|c|} \hline & \textbf{Pearson Residual} & \textbf{Standardized Residual} \\ \hline \text{Kelompok 1, Respon Positif} & 4.00 & 6.97 \\ \text{Kelompok 1, Respon Negatif} & -3.27 & -7.70 \\ \text{Kelompok 2, Respon Positif} & -3.62 & -6.30 \\ \text{Kelompok 2, Respon Negatif} & 2.96 & 6.97 \\ \hline \end{array} \]
Interpretasi :
Sel (Kelompok 1, Respon Positif) dan (Kelompok 2, Respon Negatif) menunjukkan residual Pearson positif yang besar (4.00 dan 2.96) serta standardized residual masing-masing 6.97 → Artinya, jumlah observasi aktual jauh lebih besar daripada yang diprediksi oleh model jika diasumsikan tidak ada hubungan antar variabel.
Sebaliknya, Sel (Kelompok 1, Respon Negatif) dan (Kelompok 2, Respon Positif) memiliki residual negatif signifikan (–3.27 dan –3.62) dengan standardized residual –7.70 dan –6.30 → Ini menunjukkan bahwa frekuensi aktual jauh lebih rendah dari yang diharapkan jika kategori benar-benar independen.
Karena seluruh standardized residual melebihi |3|, maka seluruh sel bisa dianggap sebagai outlier signifikan, yang memberikan sinyal kuat tentang adanya pola yang tidak acak dalam data.
Dengan demikian, dapat disimpulkan bahwa hubungan antara Kelompok dan Respon sangat tidak seimbang → Terdapat hubungan yang kuat antar kategori, dan asumsi independensi tidak terpenuhi.
Perhitungan menggunakan software r
# Data Observasi
data <- matrix(c(60, 30, 20, 90), nrow = 2, byrow = TRUE)
data
## [,1] [,2]
## [1,] 60 30
## [2,] 20 90
# Hitung nilai ekspektasi
expected <- chisq.test(data)$expected
# Pearson Residual
pearson_residual <- (data - expected) / sqrt(expected)
# Gunakan expected untuk menghitung p_row dan p_col
row_sum_exp <- rowSums(expected)
col_sum_exp <- colSums(expected)
total_sum <- sum(expected)
p_row <- row_sum_exp / total_sum
p_col <- col_sum_exp / total_sum
# Buat matriks proporsi marginal
p_row_mat <- matrix(rep(p_row, each = length(p_col)), nrow = 2)
p_col_mat <- matrix(rep(p_col, times = length(p_row)), nrow = 2, byrow = TRUE)
# Standardized Residual
standardized_residual <- (data - expected) / sqrt(expected * (1 - p_row_mat) * (1 - p_col_mat))
# Tampilkan hasil
list(Pearson_Residual = pearson_residual, Standardized_Residual = standardized_residual)
## $Pearson_Residual
## [,1] [,2]
## [1,] 4.000000 -3.265986
## [2,] -3.618136 2.954196
##
## $Standardized_Residual
## [,1] [,2]
## [1,] 6.963106 -7.698004
## [2,] -6.298367 6.963106
Analisis residual digunakan untuk mendeteksi adanya penyimpangan yang signifikan dalam tabel kontingensi, yaitu dengan melihat nilai residual dari masing-masing sel yang jauh dari nilai ekspektasi berdasarkan model independensi.
Jika nilai residual standar (standardized residual) lebih besar dari ±3, maka sel tersebut dapat dianggap sebagai indikasi adanya outlier yang signifikan yang mungkin berdampak pada kesimpulan hasil analisis.
Dari hasil yang diperoleh:
Sel (1,1) atau Kategori pertama & Sukses memiliki nilai residual standar 6.96, yang berarti jumlah sukses jauh lebih besar dari yang diharapkan.
Sel (1,2) atau Kategori pertama & Gagal memiliki residual standar -7.70, menandakan bahwa jumlah gagal jauh lebih kecil dari ekspektasi.
Sel (2,1) atau Kategori kedua & Sukses menunjukkan nilai -6.30, berarti jumlah sukses lebih kecil dari yang diprediksi.
Sel (2,2) atau Kategori kedua & Gagal memiliki nilai positif besar, yaitu 6.96, menandakan terlalu banyak kasus gagal dibandingkan harapan.
Kesimpulan sementara: Tabel kontingensi menunjukkan hubungan yang sangat kuat dan tidak acak antara kategori dan hasil, karena seluruh residual standar melebihi ±3. Hal ini menunjukkan bahwa asumsi independensi tidak terpenuhi, dan salah satu kategori jelas mendominasi dalam keberhasilan atau kegagalan.
Tabel kontingensi tiga arah merupakan alat analisis statistik lanjutan yang digunakan untuk mengevaluasi hubungan antara tiga variabel kategorik secara bersamaan. Dalam konteks ini, dua variabel utama (misalnya 𝑋 dan 𝑌) dianalisis dengan mempertimbangkan pengaruh dari variabel ketiga (𝑍) yang berfungsi sebagai variabel kontrol atau kovariat.
Analisis seperti ini sangat bermanfaat dalam berbagai bidang penelitian, khususnya ketika ada kemungkinan bahwa hubungan antara dua variabel utama dipengaruhi oleh variabel ketiga. Contoh penerapannya antara lain:
Epidemiologi: Untuk meneliti hubungan antara kebiasaan merokok (X) dan risiko kanker paru-paru (Y) dengan mempertimbangkan pengaruh usia atau jenis kelamin sebagai variabel kontrol (Z).
Ilmu sosial: Dalam mengevaluasi hubungan antara status terdakwa (X) dan hasil keputusan hukum (Y) dengan mempertimbangkan ras atau latar belakang sosial sebagai kovariat (Z).
Tabel kontingensi tiga arah biasanya dipecah menjadi dua bentuk:
Tabel kontingensi tiga arah biasanya dipecah menjadi dua bentuk:Tabel Parsial: Menunjukkan hubungan antara 𝑋 dan 𝑌dalam setiap kategori dari 𝑍. Ini memungkinkan analisis hubungan bersyarat (conditional association), di mana pengaruh 𝑍dikendalikan.
Tabel Marginal: Menyajikan hubungan antara 𝑋dan𝑌tanpa mempertimbangkan 𝑍(dengan menjumlahkan semua kategori𝑍). Analisis ini berguna untuk melihat gambaran umum, tetapi bisa menimbulkan bias interpretasi seperti pada fenomena Simpson’s Paradox.
Tujuan Penggunaan Tabel Kontingensi Tiga Arah
Mengendalikan efek kovariat yang berpotensi mengganggu hubungan antara dua variabel utama.
Menyediakan gambaran yang lebih lengkap dan akurat mengenai hubungan antar variabel.
Menghindari misinterpretasi yang bisa terjadi jika hanya melihat data secara agregat.
Manfaat Tabel Kontingensi Tiga Arah
Validitas analisis meningkat: Dengan mengendalikan variabel perancu, kesimpulan menjadi lebih valid.
Identifikasi hubungan kompleks: Mampu menangkap interaksi antar variabel yang tidak terlihat pada analisis dua arah.
Deteksi Simpson’s Paradox: Menghindari kesalahan interpretasi akibat penggabungan data yang menyembunyikan pola asli.
Keunggulan
Lebih akurat dalam analisis hubungan antar variabel.
Memungkinkan pengendalian terhadap variabel perancu.
Cocok digunakan dalam penelitian kuantitatif berskala besar.
Kekurangan
Lebih kompleks dalam perhitungan dan interpretasi.
Membutuhkan data yang cukup besar untuk menjaga kestabilan estimasi pada setiap kategori.
Tidak cocok jika variabel kontrol memiliki banyak kategori (karena dapat menimbulkan tabel yang terlalu besar dan sulit diinterpretasikan).
Tabel parsial adalah tabel yang mengelompokan X dan Y berdasarkan setiap level Z, sedangkan tabel marginal adalah tabel yang mengabaikan Z, dengan menjumlahkan data dari semua level Z.
Struktur Tabel Kontingensi Tiga Arah
Contoh kasus:
Diabetes mellitus merupakan salah satu penyakit tidak menular (PTM) yang terus meningkat secara global dan menjadi tantangan kesehatan utama abad ke-21. Berdasarkan data International Diabetes Federation (IDF), jumlah penderita diabetes di dunia diperkirakan mencapai 537 juta orang pada tahun 2021, dan angka ini diproyeksikan meningkat menjadi 643 juta pada tahun 2030 jika tidak ada upaya pencegahan yang efektif (IDF, 2021). Di Indonesia sendiri, prevalensi diabetes juga mengalami peningkatan yang signifikan, dengan data Riskesdas 2018 menunjukkan bahwa 10,9% penduduk berusia di atas 15 tahun menderita diabetes (Kemenkes RI, 2018).
Salah satu faktor risiko utama yang dikaitkan dengan diabetes adalah obesitas, yang telah terbukti meningkatkan risiko resistensi insulin dan gangguan metabolisme lainnya (Kahn et al., 2006). Namun, selain obesitas, faktor rentang usia juga memainkan peran penting dalam menentukan risiko seseorang terhadap diabetes. Beberapa studi menunjukkan bahwa risiko diabetes meningkat seiring bertambahnya usia, tetapi terdapat perbedaan pola hubungan antar kelompok usia yang belum sepenuhnya dipahami (Menke et al., 2015). Oleh karena itu, diperlukan analisis lebih lanjut yang mempertimbangkan interaksi antara usia dan obesitas dalam mempengaruhi kejadian diabetes.
Dalam analisis ini, variabel utama yang diteliti adalah hubungan antara obesitas dan diabetes. Namun, terdapat faktor lain yang dapat memengaruhi atau mengaburkan hubungan ini, salah satunya adalah usia. Usia dapat mempengaruhi baik obesitas maupun diabetes. Studi menunjukkan bahwa risiko diabetes meningkat seiring bertambahnya usia karena perubahan metabolisme, resistensi insulin, dan penurunan aktivitas fisik (Menke et al., 2015). Sementara itu, tingkat obesitas juga dapat bervariasi pada kelompok usia yang berbeda, yang pada akhirnya dapat membuat hubungan antara obesitas dan diabetes menjadi tidak sepenuhnya linear.
Data ini menunjukan hasil survei kesehatan di Indonesia yang mengkaji hubungan antara obesitas dan kejadian diabetes pada berbagai kelompok usia. Terdapat lima kelompok usia yang diamati, mulai dari balita hingga lansia. Hasil survei menunjukkan bahwa prevalensi diabetes cenderung meningkat pada individu obesitas, terutama di kelompok usia dewasa madya (50–69 tahun) dan lansia (>70 tahun), di mana jumlah penderita diabetes pada individu obesitas jauh lebih tinggi dibandingkan yang tidak obesitas. Sebaliknya, pada kelompok usia muda seperti anak-anak (5–14 tahun) dan balita, prevalensi diabetes masih rendah, baik pada individu obesitas maupun tidak. Pola ini mengindikasikan bahwa obesitas berkontribusi signifikan terhadap risiko diabetes, dan pengaruh tersebut semakin kuat seiring bertambahnya usia, sehingga penting untuk mempertimbangkan usia sebagai faktor yang memoderasi hubungan antara obesitas dan diabetes di Indonesia. Data ini digenerate menggunakan bantuan ChatGPT untuk mendukung analisis hubungan obesitas, diabetes, dan usia secara lebih komprehensif.
Tabel berikut menggambarkan hubungan antara usia (Z), Obesitas (X). dan kejadian obesitas (Y).
\[ \begin{array}{|c|c|c|c|} \hline \textbf{Usia (Tahun)} & \textbf{Obesitas} & \textbf{Diabetes: Ya} & \textbf{Diabetes: Tidak} \\ \hline < 5 & Ya & 1 & 9 \\ < 5 & Tidak & 2 & 8 \\ 5 - 14 & Ya & 5 & 9 \\ 5 - 14 & Tidak & 6 & 8 \\ 15 - 49 & Ya & 9 & 3 \\ 15 - 49 & Tidak & 4 & 6 \\ 50 - 69 & Ya & 12 & 5 \\ 50 - 69 & Tidak & 4 & 9 \\ \ge 70 & Ya & 14 & 4 \\ \ge 70 & Tidak & 5 & 9 \\ \hline \end{array} \]
Tabel frekuensi parsial menyajikan hubungan antara dua variabel kategori dalam tabel kontingensi tiga arah dengan mempertahankan satu variabel sebagai kontrol. Tabel ini membantu dalam memahami hubungan bersyarat antara variabel dalam analisis data kategori
Tabel Frekuensi parsial untuk Z (Usia) \(< 5\) tahun: \[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 1 & 9 & 10 \\ \textbf{Tidak Obesitas} & 2 & 8 & 10 \\ \hline \textbf{Total} & 3 & 17 & 20 \\ \hline \end{array} \]
Tabel Frekuensi parsial untuk Z (Usia) \(5 - 14\) tahun: \[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 5 & 9 & 14 \\ \textbf{Tidak Obesitas} & 6 & 8 & 14 \\ \hline \textbf{Total} & 11 & 17 & 28 \\ \hline \end{array} \]
Tabel Frekuensi parsial untuk Z (Usia) \(15 - 49\) tahun: \[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 9 & 3 & 12 \\ \textbf{Tidak Obesitas} & 4 & 6 & 10 \\ \hline \textbf{Total} & 13 & 9 & 22 \\ \hline \end{array} \]
Tabel Frekuensi parsial untuk Z (Usia) \(50 - 69\) tahun: \[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 12 & 5 & 17 \\ \textbf{Tidak Obesitas} & 4 & 9 & 13 \\ \hline \textbf{Total} & 16 & 14 & 30 \\ \hline \end{array} \]
Tabel Frekuensi parsial untuk Z (Usia) \(\ge 70\) tahun: \[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 14 & 4 & 18 \\ \textbf{Tidak Obesitas} & 5 & 9 & 14 \\ \hline \textbf{Total} & 19 & 13 & 32 \\ \hline \end{array} \] Implementasi dalam R :
# Data array asli
data3 <- array(c(
1, 9, # <5, Ya
2, 8, # <5, Tidak
5, 9, # 5-14, Ya
6, 8, # 5-14, Tidak
9, 3, # 15-49, Ya
4, 6, # 15-49, Tidak
12, 5, # 50-69, Ya
4, 9, # 50-69, Tidak
14, 4, # ≥70, Ya
5, 9 # ≥70, Tidak
),
dim = c(2, 2, 5),
dimnames = list(
Obesitas = c("Obesitas", "Tidak Obesitas"),
Diabetes = c("Diabetes", "Tidak Diabetes"),
Usia = c("<5", "5-14", "15-49", "50-69", "≥70")
))
# Fungsi untuk menampilkan tabel parsial dengan total
tabel_parsial_dengan_total <- function(tabel, nama_usia) {
tab <- addmargins(tabel)
cat(paste0("\nTabel Frekuensi parsian untuk Z (Usia) ", nama_usia, " tahun:\n"))
print(tab)
}
# Tampilkan semua tabel parsial
tabel_parsial_dengan_total(data3[, , "<5"], "< 5")
##
## Tabel Frekuensi parsian untuk Z (Usia) < 5 tahun:
## Diabetes
## Obesitas Diabetes Tidak Diabetes Sum
## Obesitas 1 2 3
## Tidak Obesitas 9 8 17
## Sum 10 10 20
tabel_parsial_dengan_total(data3[, , "5-14"], "5 – 14")
##
## Tabel Frekuensi parsian untuk Z (Usia) 5 – 14 tahun:
## Diabetes
## Obesitas Diabetes Tidak Diabetes Sum
## Obesitas 5 6 11
## Tidak Obesitas 9 8 17
## Sum 14 14 28
tabel_parsial_dengan_total(data3[, , "15-49"], "15 – 49")
##
## Tabel Frekuensi parsian untuk Z (Usia) 15 – 49 tahun:
## Diabetes
## Obesitas Diabetes Tidak Diabetes Sum
## Obesitas 9 4 13
## Tidak Obesitas 3 6 9
## Sum 12 10 22
tabel_parsial_dengan_total(data3[, , "50-69"], "50 – 69")
##
## Tabel Frekuensi parsian untuk Z (Usia) 50 – 69 tahun:
## Diabetes
## Obesitas Diabetes Tidak Diabetes Sum
## Obesitas 12 4 16
## Tidak Obesitas 5 9 14
## Sum 17 13 30
tabel_parsial_dengan_total(data3[, , "≥70"], "≥ 70")
##
## Tabel Frekuensi parsian untuk Z (Usia) ≥ 70 tahun:
## Diabetes
## Obesitas Diabetes Tidak Diabetes Sum
## Obesitas 14 5 19
## Tidak Obesitas 4 9 13
## Sum 18 14 32
Kesimpulan Tabel frekuensi parsial memungkinkan analisis hubungan antara dua variabel dengan mempertimbangkan variabel kontrol. Ini membantu mengidentifkasi pola yang mungkin tersembunyi dalam tabel frekuensi marginal.
Tabel Frekuensi Marginal menampilkan jumlah total observasi untuk setiap variabel dengan mengabaikan variabel lainnya dalam tabel kontingensi tiga arah. Tabel ini membantu dalam memahami distribusi kategori secara agregat tanpa mempertimbangkan hubungan antarvariabel. Tabel frekuensi marginal dihitung dengan menjumlahkan frekuensi dari tabel kontingensi tiga arah berdasarkan variabel yang tersisa.
Tabel Frekuensi Marginal untuk X (Obesitas) dan Y (Diabetes): \[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 41 & 30 & 71 \\ \textbf{Tidak Obesitas} & 21 & 40 & 61 \\ \hline \textbf{Total} & 62 & 70 & 132 \\ \hline \end{array} \]
Tabel Frekuensi Marginal untuk Z (Usia) dan Y (Diabetes):
\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline < 5 & 3 & 17 & 20 \\ 5 - 14 & 11 & 17 & 28 \\ 15 - 49 & 13 & 9 & 22 \\ 50 - 69 & 16 & 14 & 30 \\ \ge70 & 19 & 13 & 32 \\ \hline \textbf{Total} & 62 & 70 & 132 \\ \hline \end{array} \]
# Data array asli
data3 <- array(c(
1, 9, # <5, Ya
2, 8, # <5, Tidak
5, 9, # 5-14, Ya
6, 8, # 5-14, Tidak
9, 3, # 15-49, Ya
4, 6, # 15-49, Tidak
12, 5, # 50-69, Ya
4, 9, # 50-69, Tidak
14, 4, # ≥70, Ya
5, 9 # ≥70, Tidak
),
dim = c(2, 2, 5),
dimnames = list(
Obesitas = c("Obesitas", "Tidak Obesitas"),
Diabetes = c("Diabetes", "Tidak Diabetes"),
Usia = c("<5", "5-14", "15-49", "50-69", "≥70")
))
# Fungsi untuk menampilkan tabel parsial dengan total
freq_marginal_X <- apply(data3,1,sum)
freq_marginal_Z <- apply(data3,3,sum)
#TAMPILKAN HASIL
freq_marginal_X
## Obesitas Tidak Obesitas
## 62 70
freq_marginal_Z
## <5 5-14 15-49 50-69 ≥70
## 20 28 22 30 32
\[ P(Z,X,Y) = \frac{f(Z,X,Y)}{N} \]
Sebagai contoh, kita akan menghitung beberapa peluang bersama dari dari kasus yang sebelumnya.
\[ P(\text{< 5},\ \text{Obesitas},\ \text{Diabetes}) = \frac{1}{132} \approx 0.0076 \]
\[ P(\text{< 5},\ \text{Obesitas},\ \text{Tidak Diabetes}) = \frac{9}{132} \approx 0.0682 \]
\[ P(\text{< 5},\ \text{Tidak Obesitas},\ \text{Diabetes}) = \frac{2}{132} \approx 0.0152 \]
\[ P(\text{< 5},\ \text{Tidak Obesitas},\ \text{Tidak Diabetes}) = \frac{8}{132} \approx 0.0606 \]
\[ P(\text{5–14},\ \text{Obesitas},\ \text{Diabetes}) = \frac{5}{132} \approx 0.0379 \] \[ P(\text{5–14},\ \text{Obesitas},\ \text{Tidak Diabetes}) = \frac{9}{132} \approx 0.0682 \] \[ P(\text{5–14},\ \text{Tidak Obesitas},\ \text{Diabetes}) = \frac{6}{132} \approx 0.0455 \]
\[ P(\text{5–14},\ \text{Tidak Obesitas},\ \text{Tidak Diabetes}) = \frac{8}{132} \approx 0.0606 \]
Berikut tabel peluang bersama:
\[ \begin{array}{|c|c|c|c|} \hline \textbf{Usia (Tahun)} & \textbf{Obesitas} & \textbf{Diabetes: Ya} & \textbf{Diabetes: Tidak} \\ \hline < 5 & Ya & 0.0076 & 0.0682 \\ < 5 & Tidak & 0.0152 & 0.0606 \\ 5 - 14 & Ya & 0.0379 & 0.0682 \\ 5 - 14 & Tidak & 0.0455 & 0.0606 \\ 15 - 49 & Ya & 0.0682 & 0.227 \\ 15 - 49 & Tidak & 0.0303 & 0.0455 \\ 50 - 69 & Ya & 0.0909 & 0.0379 \\ 50 - 69 & Tidak & 0.0303 & 0.682 \\ \ge 70 & Ya & 0.1061 & 0.0303 \\ \ge 70 & Tidak & 0.0379 & 0.0682 \\ \hline \end{array} \] Implementasi dalam r
# Data frekuensi sesuai urutan dari gambar kamu (per baris usia)
f <- c(
1, 2, 9, 8, # < 5
5, 6, 9, 8, # 5–14
9, 4, 3, 6, # 15–49
12, 4, 5, 9, # 50–69
14, 5, 4, 9 # ≥ 70
)
# Total frekuensi
N <- sum(f) # Tetap 132
# Hitung peluang
p <- signif(f / N, 4)
# Buat urutan kategori sesuai format kamu
usia <- rep(c("< 5", "5–14", "15–49", "50–69", "≥70"), each = 4)
obesitas <- rep(c("Ya", "Tidak", "Ya", "Tidak"), 5)
diabetes <- rep(c("Ya", "Ya", "Tidak", "Tidak"), 5)
# Buat dataframe
data <- data.frame(
Usia = usia,
Obesitas = obesitas,
Diabetes = diabetes,
Frekuensi = f,
Peluang = p
)
# Lihat hasil
print(data)
## Usia Obesitas Diabetes Frekuensi Peluang
## 1 < 5 Ya Ya 1 0.007576
## 2 < 5 Tidak Ya 2 0.015150
## 3 < 5 Ya Tidak 9 0.068180
## 4 < 5 Tidak Tidak 8 0.060610
## 5 5–14 Ya Ya 5 0.037880
## 6 5–14 Tidak Ya 6 0.045450
## 7 5–14 Ya Tidak 9 0.068180
## 8 5–14 Tidak Tidak 8 0.060610
## 9 15–49 Ya Ya 9 0.068180
## 10 15–49 Tidak Ya 4 0.030300
## 11 15–49 Ya Tidak 3 0.022730
## 12 15–49 Tidak Tidak 6 0.045450
## 13 50–69 Ya Ya 12 0.090910
## 14 50–69 Tidak Ya 4 0.030300
## 15 50–69 Ya Tidak 5 0.037880
## 16 50–69 Tidak Tidak 9 0.068180
## 17 ≥70 Ya Ya 14 0.106100
## 18 ≥70 Tidak Ya 5 0.037880
## 19 ≥70 Ya Tidak 4 0.030300
## 20 ≥70 Tidak Tidak 9 0.068180
Sebagai contoh: \[ P(\text{Diabetes} = \text{Ya}) = \frac{1 + 2 + 5 + 6 + 9 + 4 + 12 + 4 + 14 + 5}{132} = \frac{62}{132} \approx 0.4697 \]
\[ P(\text{Obesitas} = \text{Tidak}) = \frac{2 + 8 + 6 + 8 + 4 + 6 + 4 + 9 + 5 + 9}{132} = \frac{61}{132} \approx 0.46212 \]
Tabel Peluang Marginal untuk X (Obesitas) dan Y (Diabetes)
\[ \begin{array}{|c|c|c|c|} \hline \textbf{X (Obesitas)}& \textbf{P(Y)} & \textbf{P (Tidak Y)} \\ \hline \textbf{Ya} & 0.31061 & 0.22727 \\ \textbf{Tidak} & 0.15909 & 0.30303 \\ \hline \end{array} \]
Tabel Peluang Marginal untuk Z (Usia) dan Y (Diabetes) \[ \begin{array}{|c|c|c|c|} \hline \textbf{X (Obesitas)}& \textbf{P(Y)} & \textbf{P (Tidak Y)} \\ \hline < 5 & 0.0227 & 0.1288 \\ 5 - 14 & 0.0833 & 0.1288 \\ 14 - 49 & 0.0985 & 0.0682 \\ 50 - 69 & 0.1212 & 0.1061 \\ \ge 70 & 0.1439 & 0.0985 \\ \hline \end{array} \]
# Data frekuensi asli (dari tabel gambar)
data3 <- array(
c(
1, 9, # <5, Ya
2, 8, # <5, Tidak
5, 9, # 5–14, Ya
6, 8, # 5–14, Tidak
9, 3, # 15–49, Ya
4, 6, # 15–49, Tidak
12, 5, # 50–69, Ya
4, 9, # 50–69, Tidak
14, 4, # ≥70, Ya
5, 9 # ≥70, Tidak
),
dim = c(2, 2, 5), # Obesitas (Ya/Tidak), Diabetes (Ya/Tidak), Usia (5 kategori)
dimnames = list(
Obesitas = c("Ya", "Tidak"),
Diabetes = c("Ya", "Tidak"),
Usia = c("<5", "5-14", "15-49", "50-69", "≥70")
)
)
# Hitung total
total <- sum(data3)
# Hitung probabilitas bersama
joint_prob <- data3 / total
# Hitung probabilitas marginal
# Untuk Obesitas
marginal_Obes <- apply(joint_prob, 1, sum)
# Untuk Diabetes
marginal_Diab <- apply(joint_prob, 2, sum)
# Untuk Usia
marginal_Usia <- apply(joint_prob, 3, sum)
# Tampilkan hasil
marginal_Obes
## Ya Tidak
## 0.469697 0.530303
marginal_Diab
## Ya Tidak
## 0.5378788 0.4621212
marginal_Usia
## <5 5-14 15-49 50-69 ≥70
## 0.1515152 0.2121212 0.1666667 0.2272727 0.2424242
Kesimpulan Tabel peluang marginal membantu memahami distribusi probabilitas dari masing-masing variabel tanpa mempertimbangkan variabel lainnya. Ini memberikan wawasan penting dalam analisis data kategori dengan tabel kontingensi tiga arah.
Peluang bersyarat didefinisikan sebagai: \[ P(Z|X,Y) = \frac{P(Z,X,Y)}{P(X,Y)} \] \[ \begin{aligned} &\text{Total populasi: } N = 132 \\ \\ &\text{Probabilitas bersama: } \\ &P(Z = \geq 70, X = \text{Ya}, Y = \text{Ya}) = \frac{14}{132} \\ &P(X = \text{Ya}, Y = \text{Ya}) = \frac{1 + 5 + 9 + 12 + 14}{132} = \frac{41}{132} \\ \\ &\text{Peluang bersyarat:} \\ &P(Z = \geq 70 \mid X = \text{Ya}, Y = \text{Ya}) = \frac{P(Z = \geq 70, X = \text{Ya}, Y = \text{Ya})}{P(X = \text{Ya}, Y = \text{Ya})} \\ &= \frac{14 / 132}{41 / 132} = \frac{14}{41} \approx 0.341 \\ \end{aligned} \]
Implementasi dalam R
# Frekuensi data
data3 <- array(
c(
1, 9, # <5, Ya
2, 8, # <5, Tidak
5, 9, # 5–14, Ya
6, 8, # 5–14, Tidak
9, 3, # 15–49, Ya
4, 6, # 15–49, Tidak
12, 5, # 50–69, Ya
4, 6, # 50–69, Tidak
14, 4, # ≥70, Ya
5, 9 # ≥70, Tidak
),
dim = c(2, 2, 5),
dimnames = list(
Obesitas = c("Ya", "Tidak"),
Diabetes = c("Ya", "Tidak"),
Usia = c("<5", "5-14", "15-49", "50-69", "≥70")
)
)
# Hitung total populasi
total <- sum(data3)
# Hitung probabilitas bersama
joint_prob <- data3 / total
# P(X = Ya, Y = Ya)
P_XY_YaYa <- sum(joint_prob["Ya", "Ya", ])
# P(Z = ≥70, X = Ya, Y = Ya)
P_Z70_XYa_YYa <- joint_prob["Ya", "Ya", "≥70"]
# Peluang bersyarat
P_Z70_given_XYa_YYa <- P_Z70_XYa_YYa / P_XY_YaYa
# Tampilkan hasil
P_Z70_given_XYa_YYa
## [1] 0.3414634
Dalam studi epidemiologi, analisis hubungan antara dua variabel biner seperti paparan (misalnya obesitas) dan hasil (misalnya diabetes) seringkali dilakukan menggunakan ukuran asosiasi seperti Risk Difference (RD), Relative Risk (RR), dan Odds Ratio (OR). Namun, hubungan ini bisa dipengaruhi oleh variabel ketiga (misalnya usia), sehingga diperlukan analisis menggunakan tabel kontingensi 3 arah.
Menurut Agresti (2002) dalam bukunya Categorical Data Analysis, “analisis kontingensi 3 arah sangat penting untuk mengevaluasi apakah asosiasi antara dua variabel tetap stabil setelah mengendalikan variabel ketiga”. Pendekatan ini memungkinkan analisis asosiasi bersyarat (conditional association).
Definisi Ukuran Asosiasi
Risk Difference mengukur selisih risiko antara kelompok yang terpapar dan tidak terpapar.
\[ RD = P(Y=1|X=1, Z=z) - P(Y=1|X=0, Z=z) \]
dimana: - \(X\): variabel paparan (misal: obesitas), - \(Y\): variabel outcome (misal: diabetes), - \(Z\): variabel ketiga (misal: usia)
Relative Risk membandingkan proporsi kejadian outcome antar dua kelompok paparan.
\[ RR = \frac{P(Y=1|X=1, Z=z)}{P(Y=1|X=0, Z=z)} \]
Odds Ratio membandingkan odds (peluang relatif) terjadinya outcome antar dua kelompok.
\[ OR = \frac{ \frac{P(Y=1|X=1, Z=z)}{P(Y=0|X=1, Z=z)} }{ \frac{P(Y=1|X=0, Z=z)}{P(Y=0|X=0, Z=z)} } \]
Atau dalam bentuk tabel 2x2 bersyarat:
\[ OR = \frac{a \cdot d}{b \cdot c} \]
dimana: - \(a = \text{frekuensi}(X=1, Y=1, Z=z)\) - \(b = \text{frekuensi}(X=1, Y=0, Z=z)\) - \(c = \text{frekuensi}(X=0, Y=1, Z=z)\) - \(d = \text{frekuensi}(X=0, Y=0, Z=z)\)
Contoh perhitungan
Kelompok Usia: \(\le 5\) Tahun
\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 1 & 9 & 10 \\ \textbf{Tidak Obesitas} & 2 & 8 & 10 \\ \hline \textbf{Total} & 3 & 17 & 20 \\ \hline \end{array} \] Perhitungan Manual:
Risk Difference:
\[ RD = \frac{1}{1 + 9} - \frac{2}{2 + 8} = 0.1 - 0.2 = -0.1 \]
Relative Risk:
\[ RR = \frac{0.1}{0.2} = 0.5 \]
Odds Ratio:
\[ OR = \frac{1 \times 8}{9 \times 2} = \frac{8}{18} \approx 0.444 \]
Perhitungan di R:
a <- 1; b <- 9; c <- 2; d <- 8
risk_exp <- a / (a + b)
risk_unexp <- c / (c + d)
RD <- risk_exp - risk_unexp
RR <- risk_exp / risk_unexp
OR <- (a * d) / (b * c)
list(RD = RD, RR = RR, OR = OR)
## $RD
## [1] -0.1
##
## $RR
## [1] 0.5
##
## $OR
## [1] 0.4444444
Interpretasi:
Pada kelompok usia \(\le5\) tahun,
obesitas tampaknya berhubungan negatif terhadap diabetes (RR = 0.5, OR ≈
0.444).
Kelompok Usia: \(5 - 14\) Tahun
\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 5 & 9 & 14 \\ \textbf{Tidak Obesitas} & 6 & 8 & 14 \\ \hline \textbf{Total} & 11 & 17 & 28 \\ \hline \end{array} \]
Perhitungan Manual:
Risk Difference:
\[ RD = \frac{5}{14} - \frac{6}{14} = 0.357 - 0.429 = -0.072 \]
Relative Risk:
\[ RR = \frac{0.357}{0.429} \approx 0.833 \]
Odds Ratio:
\[ OR = \frac{5 \times 8}{9 \times 6} = \frac{40}{54} \approx 0.741 \]
Perhitungan di R:
a <- 5; b <- 9; c <- 6; d <- 8
risk_exp <- a / (a + b)
risk_unexp <- c / (c + d)
RD <- risk_exp - risk_unexp
RR <- risk_exp / risk_unexp
OR <- (a * d) / (b * c)
list(RD = RD, RR = RR, OR = OR)
## $RD
## [1] -0.07142857
##
## $RR
## [1] 0.8333333
##
## $OR
## [1] 0.7407407
Interpretasi:
Risiko diabetes sedikit lebih rendah pada anak obesitas (RR ≈ 0.833, OR
≈ 0.741).
Kelompok Usia: \(15 - 49\) Tahun
\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 9 & 3 & 12 \\ \textbf{Tidak Obesitas} & 4 & 6 & 10 \\ \hline \textbf{Total} & 13 & 9 & 22 \\ \hline \end{array} \] Perhitungan Manual:
Risk Difference:
\[ RD = \frac{9}{12} - \frac{4}{10} = 0.75 - 0.4 = 0.35 \]
Relative Risk:
\[ RR = \frac{0.75}{0.4} = 1.875 \]
Odds Ratio:
\[ OR = \frac{9 \times 6}{3 \times 4} = 4.5 \]
Perhitungan di R:
a <- 9; b <- 3; c <- 4; d <- 6
risk_exp <- a / (a + b)
risk_unexp <- c / (c + d)
RD <- risk_exp - risk_unexp
RR <- risk_exp / risk_unexp
OR <- (a * d) / (b * c)
list(RD = RD, RR = RR, OR = OR)
## $RD
## [1] 0.35
##
## $RR
## [1] 1.875
##
## $OR
## [1] 4.5
Interpretasi:
Risiko diabetes pada obesitas cukup tinggi di usia 15-49 (RR ≈ 1.875, OR
= 4.5).
Kelompok Usia: \(50 - 69\) Tahun \[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 12 & 5 & 17 \\ \textbf{Tidak Obesitas} & 4 & 9 & 13 \\ \hline \textbf{Total} & 16 & 14 & 30 \\ \hline \end{array} \] Perhitungan Manual:
Risk Difference:
\[ RD = \frac{12}{17} - \frac{4}{13} \approx 0.706 - 0.308 = 0.398 \]
Relative Risk:
\[ RR = \frac{0.706}{0.308} \approx 2.29 \]
Odds Ratio:
\[ OR = \frac{12 \times 9}{5 \times 4} = 5.4 \]
Perhitungan di R:
a <- 12; b <- 5; c <- 4; d <- 9
risk_exp <- a / (a + b)
risk_unexp <- c / (c + d)
RD <- risk_exp - risk_unexp
RR <- risk_exp / risk_unexp
OR <- (a * d) / (b * c)
list(RD = RD, RR = RR, OR = OR)
## $RD
## [1] 0.39819
##
## $RR
## [1] 2.294118
##
## $OR
## [1] 5.4
Interpretasi:
Risiko diabetes sangat tinggi pada obesitas (RR ≈ 2.29, OR = 5.4) untuk
usia 50-69 tahun.
Kelompok Usia: \(\ge 70\) Tahun
\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 14 & 4 & 18 \\ \textbf{Tidak Obesitas} & 5 & 9 & 14 \\ \hline \textbf{Total} & 19 & 13 & 32 \\ \hline \end{array} \] Perhitungan Manual:
Risk Difference:
\[ RD = \frac{14}{18} - \frac{5}{14} \approx 0.778 - 0.357 = 0.421 \]
Relative Risk:
\[ RR = \frac{0.778}{0.357} \approx 2.18 \]
Odds Ratio:
\[ OR = \frac{14 \times 9}{4 \times 5} = 6.3 \]
Perhitungan di R:
a <- 14; b <- 4; c <- 5; d <- 9
risk_exp <- a / (a + b)
risk_unexp <- c / (c + d)
RD <- risk_exp - risk_unexp
RR <- risk_exp / risk_unexp
OR <- (a * d) / (b * c)
list(RD = RD, RR = RR, OR = OR)
## $RD
## [1] 0.4206349
##
## $RR
## [1] 2.177778
##
## $OR
## [1] 6.3
Interpretasi:
Pada usia $ $ tahun, obesitas meningkatkan risiko diabetes signifikan
(RR ≈ 2.18, OR = 6.3).
Kesimpulan
Dari seluruh kelompok usia, obesitas berhubungan kuat dengan peningkatan risiko diabetes di usia 15 tahun ke atas, terutama signifikan pada usia 50 tahun ke atas. Hasil perhitungan manual dan verifikasi R menunjukkan konsistensi.
Perhitungan Tabel Marginal (Keseluruhan Usia)
\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Diabetes} & \textbf{Tidak Diabetes} & \textbf{Total} \\ \hline \textbf{Obesitas} & 41 & 30 & 71 \\ \textbf{Tidak Obesitas} & 21 & 40 & 61 \\ \hline \textbf{Total} & 62 & 70 & 132 \\ \hline \end{array} \]
Perhitungan Manual:
Risk Difference:
\[ RD = \frac{41}{41 + 30} - \frac{21}{21 + 40} = \frac{41}{71} - \frac{21}{61} \approx 0.577 - 0.344 = 0.233 \]
Relative Risk:
\[ RR = \frac{0.577}{0.344} \approx 1.678 \]
Odds Ratio:
\[ OR = \frac{41 \times 40}{30 \times 21} = \frac{1640}{630} \approx 2.603 \]
Perhitungan di R:
a <- 41; b <- 30; c <- 21; d <- 40
risk_exp <- a / (a + b)
risk_unexp <- c / (c + d)
RD <- risk_exp - risk_unexp
RR <- risk_exp / risk_unexp
OR <- (a * d) / (b * c)
list(RD = RD, RR = RR, OR = OR)
## $RD
## [1] 0.2332025
##
## $RR
## [1] 1.677398
##
## $OR
## [1] 2.603175
Interpretasi: Secara keseluruhan tanpa memperhatikan kelompok usia, obesitas meningkatkan risiko diabetes secara signifikan (RR ≈ 1.678, OR ≈ 2.603, RD ≈ 0.233).
Kesimpulan
Dari seluruh kelompok usia, obesitas berhubungan kuat dengan peningkatan risiko diabetes di usia 15 tahun ke atas, terutama signifikan pada usia 50 tahun ke atas. Hasil perhitungan manual dan verifikasi R menunjukkan konsistensi.
Secara keseluruhan (tabel marginal), obesitas meningkatkan risiko diabetes dengan RR ≈ 1.678 dan OR ≈ 2.603.
Conditional independence (kemandirian bersyarat) dalam tabel kontingensi terjadi ketika dua variabel menjadi independen setelah dikendalikan oleh variabel ketiga.
Secara matematis, dua variabel \(X\) dan \(Y\) dikatakan independen secara kondisional terhadap variabel \(Z\) jika:
\[ P(X,Y|Z) = P(X|Z) \cdot P(Y|Z) \]
Dengan kata lain, probabilitas gabungan antara \(X\) dan \(Y\) setelah dikendalikan oleh \(Z\), dapat diuraikan menjadi hasil perkalian probabilitas masing-masing secara bersyarat terhadap \(Z\).
Bentuk frekuensi dalam konteks data berbasis frekuensi (tabel kontingensi 3 arah), persamaan tersebut dapat dinyatakan sebagai:
\[ \frac{n_{ijk}}{n_{k++}} = \left( \frac{n_{i+k}}{n_{k++}} \right) \times \left( \frac{n_{+jk}}{n_{k++}} \right) \]
Keterangan: - \(n_{ijk}\): frekuensi kasus pada kategori ke-\(i\) untuk \(X\), kategori ke-\(j\) untuk \(Y\), dan kategori ke-\(k\) untuk \(Z\), - \(n_{i+k}\): jumlah kasus pada \(X = i\) dan \(Z = k\) (menjumlahkan semua \(j\)), - \(n_{+jk}\): jumlah kasus pada \(Y = j\) dan \(Z = k\) (menjumlahkan semua \(i\)), - \(n_{k++}\): total kasus untuk kategori \(Z = k\) (menjumlahkan semua \(i\) dan \(j\)).
Jika persamaan di atas berlaku untuk semua kombinasi \(i, j, k\), maka kita katakan bahwa \(X \perp Y \, | \, Z\) (X dan Y kondisional independen terhadap Z).
Contoh Kasus Dengan data yang sama seperti sebelumnya \[ \begin{array}{|c|c|c|c|} \hline \textbf{Usia (Tahun)} & \textbf{Obesitas} & \textbf{Diabetes: Ya} & \textbf{Diabetes: Tidak} \\ \hline < 5 & Ya & 1 & 9 \\ < 5 & Tidak & 2 & 8 \\ 5 - 14 & Ya & 5 & 9 \\ 5 - 14 & Tidak & 6 & 8 \\ 15 - 49 & Ya & 9 & 3 \\ 15 - 49 & Tidak & 4 & 6 \\ 50 - 69 & Ya & 12 & 5 \\ 50 - 69 & Tidak & 4 & 9 \\ \ge 70 & Ya & 14 & 4 \\ \ge 70 & Tidak & 5 & 9 \\ \hline \end{array} \]
library(knitr)
## Warning: package 'knitr' was built under R version 4.3.3
library(kableExtra)
## Warning: package 'kableExtra' was built under R version 4.3.3
# Buat array data kontingensi 3 arah: Usia x Obesitas x Diabetes
# Format: [Obesitas, Diabetes, Usia]
data_obesitas_diabetes <- array(
c(
1, 9, # < 5: Ya
2, 8, # < 5: Tidak
5, 9, # 5-14: Ya
6, 8, # 5-14: Tidak
9, 3, # 15-49: Ya
4, 6, # 15-49: Tidak
12, 5, # 50-69: Ya
4, 9, # 50-69: Tidak
14, 4, # >=70: Ya
5, 9 # >=70: Tidak
),
dim = c(2, 2, 5),
dimnames = list(
"Obesitas" = c("Ya", "Tidak"),
"Diabetes" = c("Ya", "Tidak"),
"Usia" = c("< 5", "5–14", "15–49", "50–69", "≥70")
)
)
Pengujian Contional Independence Untuk menguji apakah obesitas dan diabetes independen secara kondisional terhadap usia, kita gunakan uji chi-square dalam setiap kelompok usia.
Perhitungan manual untuk setiap kelompok usia (Z) Kita akan menguji apakah:
\[ P(X, Y | Z) \approx P(X | Z) \cdot P(Y | Z) \]
dengan:
Contoh untuk Usia 50 – 69 tahun
Dari data diketahui:
Jumlah total untuk usia 50 – 69:
\[ n = 12 + 5 + 4 + 9 = 30 \]
Hitung:
\[ P(Diabetes\ |\ Obesitas = \text{Ya},\ Usia = 50\!-\!69) = \frac{12}{12 + 5} = \frac{12}{17} \approx 0.7059 \]
\[ P(Diabetes\ |\ Obesitas = \text{Tidak},\ Usia = 50\!-\!69) = \frac{4}{4 + 9} = \frac{4}{13} \approx 0.3077 \]
\[ P(Diabetes\ |\ Usia = 50\!-\!69) = \frac{12 + 4}{30} = \frac{16}{30} \approx 0.5333 \]
Interpretasi :
Karena nilai \(P(Diabetes\ |\ Obesitas =
\text{Ya},\ Usia = 50\!-\!69) \approx 0.71\) dan
\(P(Diabetes\ |\ Obesitas = \text{Tidak},\
Usia = 50\!-\!69) \approx 0.31\)
berbeda jauh dari \(P(Diabetes\ |\ Usia =
50\!-\!69) \approx 0.53\), maka: \[
P(Diabetes\ |\ Obesitas,\ Usia) \not\approx P(Diabetes\ |\ Usia)
\] Sehingga Obesitas dan Diabetes tidak independen secara
kondisional terhadap Usia = 50 – 69.
Perhitungan dengan R
# Uji Chi-Square untuk masing-masing strata usia dari data obesitas-diabetes
chisq_list <- list()
for (z in dimnames(data_obesitas_diabetes)$Usia) {
tbl <- data_obesitas_diabetes[, , z]
test_result <- chisq.test(tbl)
chisq_list[[z]] <- test_result
cat("\n=== Hasil Uji Chi-Square untuk Usia:", z, "===\n")
print(test_result)
}
## Warning in chisq.test(tbl): Chi-squared approximation may be incorrect
##
## === Hasil Uji Chi-Square untuk Usia: < 5 ===
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tbl
## X-squared = 0, df = 1, p-value = 1
##
##
## === Hasil Uji Chi-Square untuk Usia: 5–14 ===
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tbl
## X-squared = 0, df = 1, p-value = 1
## Warning in chisq.test(tbl): Chi-squared approximation may be incorrect
##
## === Hasil Uji Chi-Square untuk Usia: 15–49 ===
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tbl
## X-squared = 1.5058, df = 1, p-value = 0.2198
##
##
## === Hasil Uji Chi-Square untuk Usia: 50–69 ===
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tbl
## X-squared = 3.2294, df = 1, p-value = 0.07233
##
##
## === Hasil Uji Chi-Square untuk Usia: ≥70 ===
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tbl
## X-squared = 4.1643, df = 1, p-value = 0.04129
\[ \begin{array}{|c|c|l|} \hline \textbf{Usia} & \textbf{p-value} & \textbf{Kesimpulan} \\ \hline < 5 & 1 & \text{Tidak ada asosiasi (independen)} \\ 5{-}14 & 1 & \text{Tidak ada asosiasi (independen)} \\ 15{-}49 & 0.2198 & \text{Tidak signifikan, asosiasi tidak terbukti} \\ 50{-}69 & 0.07233 & \text{Mendekati signifikan, indikasi asosiasi lemah} \\ \geq 70 & 0.04129 & \text{Signifikan, terdapat asosiasi obesitas--diabetes} \\ \hline \end{array} \]
Dalam analisis data kategorik, sering kali kita dihadapkan pada situasi di mana terdapat dua variabel utama yang ingin dianalisis hubungannya (misalnya variabel X dan Y), tetapi hubungan tersebut dapat dipengaruhi oleh variabel ketiga yang bersifat sebagai confounding atau variabel kontrol (Z). Untuk menangani situasi seperti ini, digunakanlah tabel kontingensi tiga arah, yang memungkinkan kita mengevaluasi asosiasi antara dua variabel sembari mengontrol efek dari variabel ketiga.
Sebagai contoh, kita ingin mengetahui apakah terdapat hubungan antara kebiasaan merokok (X) dengan kejadian kanker paru-paru (Y), dengan mempertimbangkan usia sebagai variabel pengganggu (Z). Pendekatan yang digunakan adalah dengan membagi data ke dalam tabel parsial (yaitu tabel 2 × 2) untuk setiap strata dari variabel Z (misalnya, kelompok usia tertentu), serta membuat sebuah tabel marginal yang mengabaikan variabel Z.
Untuk menilai kekuatan asosiasi dalam setiap tabel parsial, digunakan ukuran asosiasi yang umum yaitu Odds Ratio (OR). Odds Ratio didefinisikan sebagai:
\[ OR = \frac{a \cdot d}{b \cdot c} \]
dengan:
Nilai OR memberikan informasi mengenai kekuatan dan arah hubungan antara dua variabel. OR > 1 menunjukkan adanya asosiasi positif, OR < 1 menunjukkan asosiasi negatif, dan OR = 1 menunjukkan tidak ada asosiasi.
Ketika kita memperoleh beberapa tabel parsial (misalnya berdasarkan kelompok usia), dan jika nilai Odds Ratio dari masing-masing tabel relatif konsisten, maka kita dapat menghitung Odds Ratio gabungan menggunakan estimasi Mantel-Haenszel (MH). Estimasi ini memberikan ringkasan yang memperhitungkan variasi antar strata.
Rumus Mantel-Haenszel Odds Ratio adalah sebagai berikut:
\[ \text{OR}_{MH} = \frac{\sum_i \frac{a_i d_i}{n_i}}{\sum_i \frac{b_i c_i}{n_i}} \]
dengan \(a_i, b_i, c_i, d_i\) merupakan isi dari tabel parsial ke-\(i\), dan \(n_i\) adalah total observasi dalam tabel tersebut.
Metode Mantel-Haenszel juga memungkinkan kita melakukan pengujian hipotesis nol bahwa tidak ada asosiasi antara X dan Y dalam semua strata Z, melalui statistik Chi-Square Mantel-Haenszel.
Dalam analisis statistik terhadap data kategorik, independensi bersyarat merupakan konsep penting yang memungkinkan kita memahami hubungan antara dua variabel utama dengan mempertimbangkan adanya variabel ketiga yang berpotensi memengaruhi hubungan tersebut. Tabel kontingensi tiga arah digunakan untuk menyelidiki hubungan antara dua variabel (misalnya, X dan Y), dalam konteks strata dari variabel kontrol ketiga (Z).
Definisi Konseptual dan Matematis Independensi bersyarat terjadi ketika dua variabel, misalnya \(X\) dan \(Y\), tidak saling bergantung satu sama lain pada setiap tingkat (strata) dari variabel ketiga \(Z\). Dalam konteks tabel kontingensi, ini berarti bahwa setelah kita mengendalikan (mengontrol) pengaruh dari \(Z\), hubungan antara \(X\) dan \(Y\) menjadi netral — tidak ada asosiasi yang terdeteksi.
Secara matematis, independensi bersyarat dituliskan sebagai:
\[ OR(X, Y \mid Z) = 1 \]
Artinya, rasio odds (odds ratio) antara \(X\) dan \(Y\) adalah 1 pada setiap strata dari variabel \(Z\), yang menunjukkan tidak adanya hubungan atau ketergantungan antara keduanya dalam konteks strata tersebut.
Walaupun \(X\) dan \(Y\) mungkin independen dalam setiap strata \(Z\), bukan berarti mereka juga independen secara keseluruhan (pada data marginal, yaitu jika \(Z\) diabaikan). Hal ini dikenal sebagai Paradox Simpson, yaitu fenomena di mana arah atau kekuatan hubungan antara dua variabel berubah jika dilihat secara keseluruhan dibandingkan dengan saat dibedakan berdasarkan strata.
Pengujian independensi bersyarat biasanya dilakukan menggunakan pendekatan statistik seperti uji Cochran-Mantel-Haenszel (CMH). Uji ini dirancang untuk menguji apakah ada asosiasi antara dua variabel biner dalam beberapa strata dari variabel ketiga.
Asumsi independensi bersyarat juga sering digunakan dalam bidang lain, seperti:
Contoh dalam aplikasi nyata misalnya, dalam studi tentang hubungan antara obesitas (X) dan diabetes (Y), kita mungkin ingin mengetahui apakah hubungan tersebut tetap signifikan setelah memperhitungkan usia (Z) sebagai variabel kontrol. Jika hasil menunjukkan bahwa OR pada setiap kelompok usia mendekati 1, maka kita dapat menyimpulkan bahwa obesitas dan diabetes bersifat independen bersyarat terhadap usia.
Namun, apabila OR bervariasi secara signifikan antar strata usia, maka asumsi independensi bersyarat tidak berlaku, dan kita harus mempertimbangkan interaksi yang lebih kompleks.
Uji Cochran-Mantel-Haenszel digunakan untuk menguji hubungan antara dua variabel kategorik sambil mengendalikan efek dari variabel perancu.
Manfaat Uji CMH - Mengontrol efek perancu dalam tabel kontingensi berlapis - Menguji hipotesis independensi bersyarat - Mengatasi confounding dalam studi observasional atau eksperimen sosial
CMH membagi data ke dalam tabel-tabel \(2 \times 2\) sesuai dengan strata variabel kontrol (\(Z\)), dan menguji hubungan antara \(X\) dan \(Y\) setelah mengendalikan \(Z\).
Hipotesis:
Statistik Uji CMH:
\[ CMH = \frac{\left[ \sum_k (n_{11k} - \mu_{11k}) \right]^2}{\sum_k \text{var}(n_{11k})} \]
Keterangan:
\[ \mu_{11k} = E(n_{11k}) = \frac{n_{1\cdot k} \cdot n_{\cdot 1k}}{n_{\cdot \cdot k}} \]
\[ \text{var}(n_{11k}) = \frac{n_{1\cdot k} \cdot n_{2\cdot k} \cdot n_{\cdot 1k} \cdot n_{\cdot 2k}}{n_{\cdot\cdot k}^2(n_{\cdot\cdot k}-1)} \]
Statistik CMH mengikuti distribusi Chi-Square dengan derajat kebebasan 1.
Tolak \(H_0\) jika:
\[ CMH > \chi^2_{(1)} \quad \text{atau p-value} < \alpha \]
Contoh Kasus Kita ingin menguji apakah obesitas dan diabetes bersifat independen secara kondisional setelah mengendalikan usia menggunakan Uji Cochran-Mantel-Haenszel (CMH).
Hipotesis :
H₀ (Hipotesis Nol): Obesitas dan diabetes independen secara kondisional setelah mengendalikan variabel usia.
H₁ (Hipotesis Alternatif): Obesitas dan diabetes tidak independen secara kondisional setelah mengendalikan variabel usia.
Rumus nilai ekspektasi:
\[ \mu_{11k} = \frac{n_{1\cdot k} \cdot n_{\cdot 1k}}{n_{\cdot\cdot k}} \]
\[ \begin{aligned} \mu_{11(1)} &= \frac{10 \times 3}{20} = 1.5 \\ \mu_{11(2)} &= \frac{14 \times 11}{28} = 5.5 \\ \mu_{11(3)} &= \frac{12 \times 13}{22} = 7.09 \\ \mu_{11(4)} &= \frac{17 \times 16}{30} = 9.07 \\ \mu_{11(5)} &= \frac{18 \times 19}{32} = 10.69 \\ \end{aligned} \]
Rumus varians:
\[ \text{Var}(n_{11k}) = \frac{n_{1\cdot k} \cdot n_{2\cdot k} \cdot n_{\cdot 1k} \cdot n_{\cdot 2k}}{n_{\cdot\cdot k}^2 (n_{\cdot\cdot k} - 1)} \]
\[ \begin{aligned} \text{Var}(n_{11(1)}) &= \frac{10 \times 10 \times 3 \times 17}{20^2 \times 19} = 0.67105 \\ \text{Var}(n_{11(2)}) &= \frac{14 \times 14 \times 11 \times 17}{28^2 \times 27} = 1.73148 \\ \text{Var}(n_{11(3)}) &= \frac{12 \times 10 \times 13 \times 9}{22^2 \times 21} = 1.38135 \\ \text{Var}(n_{11(4)}) &= \frac{17 \times 13 \times 16 \times 14}{30^2 \times 29} = 1.8967 \\ \text{Var}(n_{11(5)}) &= \frac{18 \times 14 \times 19 \times 13}{32^2 \times 31} = 1.96081 \\ \end{aligned} \]
Menghitung Statistik Uji CMH
\[ X^2_{CMH} = \frac{\left( \sum_k (n_{11k} - \mu_{11k}) \right)^2}{\sum_k \text{Var}(n_{11k})} \]
\[ \begin{aligned} X^2_{CMH} &= \frac{51.1929}{0,67105 + 1,73148 + 1,38135 + 1,8967 + 1,96081} \\ &= \frac{51.1929}{7.6414} \\ &= 6.69942 \end{aligned} \] Statistik uji ini dibandingkan dengan distribusi \(\chi^2\) dengan 1 derajat kebebasan. Jika nilai \(X^2_{CMH}\) lebih besar dari nilai kritis \(\chi^2_{1, \alpha}\), maka kita tolak hipotesis nol dan menyimpulkan bahwa ada hubungan antara obesitas dan diabetes setelah mengontrol variabel usia.
Dengan \(\alpha = 0.05\), nilai kritis \(\chi^2_{1,0.05} = 3.841\). Karena \(6.69942 > 3.841\), maka hipotesis nol ditolak, yang berarti ada hubungan signifikan antara obesitas dan diabetes setelah mengendalikan variabel usia. Ini menunjukkan bahwa obesitas dan diabetes tidak bersifat conditional independence setelah usia dikendalikan.
Perhitungan CMH secara manual memerlukan substitusi data ke dalam rumus di atas. Untuk memastikan hasil, kita akan menggunakan R.
Implementasi dalam R
# Data dalam bentuk array 3D (Usia, Obesitas, Diabetes)
data <- array(
c(
1, 2, # Usia <5, Obesitas=Ya (Diabetes=Ya, Tidak)
9, 8, # Usia <5, Obesitas=Tidak
5, 6, # Usia 5–14, Obesitas=Ya
9, 8, # Usia 5–14, Obesitas=Tidak
9, 4, # Usia 15–49, Obesitas=Ya
3, 6, # Usia 15–49, Obesitas=Tidak
12, 4, # Usia 50–69, Obesitas=Ya
5, 9, # Usia 50–69, Obesitas=Tidak
14, 5, # Usia ≥70, Obesitas=Ya
4, 9 # Usia ≥70, Obesitas=Tidak
),
dim = c(2, 2, 5),
dimnames = list(
Obesitas = c("Ya", "Tidak"),
Diabetes = c("Ya", "Tidak"),
Usia = c("<5", "5–14", "15–49", "50–69", "≥70")
)
)
data
## , , Usia = <5
##
## Diabetes
## Obesitas Ya Tidak
## Ya 1 9
## Tidak 2 8
##
## , , Usia = 5–14
##
## Diabetes
## Obesitas Ya Tidak
## Ya 5 9
## Tidak 6 8
##
## , , Usia = 15–49
##
## Diabetes
## Obesitas Ya Tidak
## Ya 9 3
## Tidak 4 6
##
## , , Usia = 50–69
##
## Diabetes
## Obesitas Ya Tidak
## Ya 12 5
## Tidak 4 9
##
## , , Usia = ≥70
##
## Diabetes
## Obesitas Ya Tidak
## Ya 14 4
## Tidak 5 9
# Memuat paket yang dibutuhkan
library(stats)
# Menjalankan uji Cochran-Mantel-Haenszel
cmh_test <- mantelhaen.test(data)
# Menampilkan hasil
print(cmh_test)
##
## Mantel-Haenszel chi-squared test with continuity correction
##
## data: data
## Mantel-Haenszel X-squared = 5.7958, df = 1, p-value = 0.01606
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 1.232444 5.208117
## sample estimates:
## common odds ratio
## 2.533518
Kriteria uji:
Jika p-value < 0.05, maka kita menolak H₀, yang berarti obesitas dan diabetes tidak independen secara kondisional setelah mengendalikan usia. Jika p-value > 0.05, maka kita gagal menolak H₀, sehingga obesitas dan diabetes dapat dianggap independen secara kondisional setelah mengendalikan usia.
Keputusan dan Kesimpulan
Karena p-value = 0.01606449 < 0.05, maka hipotesis nol ditolak, yang berarti ada hubungan signifikan antara obesitas dan diabetes setelah mengendalikan variabel usia. Dengan kata lain, obesitas dan diabetes tidak bersifat conditional independence setelah usia dikendalikan.
Odds Ratio (OR) adalah ukuran asosiasi dalam tabel kontingensi 2x2 yang menunjukkan seberapa besar peluang suatu kejadian terjadi pada satu kelompok dibandingkan dengan kelompok lain. Namun, ketika data dikategorikan menjadi beberapa strata berdasarkan variabel perancu (seperti usia), kita tidak cukup hanya menghitung odds ratio dari satu tabel. Untuk itu, kita gunakan Odds Ratio Bersama atau Common Odds Ratio menggunakan pendekatan Mantel-Haenszel (MH).
Odds ratio bersama (ditaksir oleh Mantel-Haenszel) dirumuskan sebagai berikut:
\[ \hat{\theta}_{MH} = \frac{\sum_{k=1}^{K} \left( \frac{n_{11k} n_{22k}}{n_{..k}} \right)}{\sum_{k=1}^{K} \left( \frac{n_{12k} n_{21k}}{n_{..k}} \right)} \] - \(n_{11k}\) = frekuensi sel baris 1 kolom 1 pada tabel parsial ke-\(k\) (misalnya: Obesitas = “Ya”, Diabetes = “Ya”) - \(n_{12k}\) = frekuensi sel baris 1 kolom 2 pada tabel parsial ke-\(k\) (Obesitas = “Ya”, Diabetes = “Tidak”) - \(n_{21k}\) = frekuensi sel baris 2 kolom 1 pada tabel parsial ke-\(k\) (Obesitas = “Tidak”, Diabetes = “Ya”) - \(n_{22k}\) = frekuensi sel baris 2 kolom 2 pada tabel parsial ke-\(k\) (Obesitas = “Tidak”, Diabetes = “Tidak”) - \(n_{\cdot\cdot k}\) = total pengamatan dalam tabel parsial ke-\(k\)
Standar Error log Odds Ratio Bersama:
\[ \hat{\sigma}^{2}[\log(\hat{\theta}_{MH})] = \frac{\sum (n_{11k} + n_{12k})(n_{11k}n_{22k})/n_{..k}^{2}}{2(\sum n_{11k}n_{12k}/n_{..k})^{2}} + \frac{\sum(n_{11k} + n_{22k})(n_{11k} + n_{12k}) + (n_{12k} + n_{21k})(n_{11k} + n_{22k})]/n_{..k}^{2}}{2(\sum n_{11k}n_{12k}/n_{..k})^{2}} + \frac{\sum(n_{12k} + n_{21k})(n_{12k}n_{21k})/n_{..k}^{2}}{2(\sum n_{12k}n_{21k}/n_{..k})^{2}} \]
Interval Kepercayaan log Odds Ratio Bersama:
\[ \log(\hat{\theta}_{MH}) \pm Z_{\alpha} \hat{\sigma}[\log(\hat{\theta}_{MH})] \] Implementasikan ke dalam kasus sebelumnya mengenai obesitas dan diabetes
\[ \begin{align*} \text{Pembilang} &= \frac{(1)(8)}{20} + \frac{(5)(8)}{28} + \frac{(9)(6)}{22} + \frac{(12)(9)}{30} + \frac{(14)(9)}{32} \\ &= 0.4 + 1.43 + 2.45 + 3.6 + 3.94 = 11.82 \\ \\ \text{Penyebut} &= \frac{(9)(2)}{20} + \frac{(9)(6)}{28} + \frac{(3)(4)}{22} + \frac{(5)(4)}{30} + \frac{(4)(5)}{32} \\ &= 0.9 + 1.93 + 0.55 + 0.67 + 0.63 = 4.68 \\ \\ \hat{\theta}_{MH} &= \frac{11.82}{4.68} = 2.53 \end{align*} \]
Misalkan simpangan baku dari \(\log(\hat{\theta}_{MH})\) diketahui sebesar 0.05:
\[ SE(\log \hat{\theta}_{MH}) = 0.05 \]
\[ \begin{align*} \log(2.53) &\pm 1.96 \times 0.05 \\ 0.927 &\pm 0.098 \\ \Rightarrow \text{CI log} &= [0.829,\ 1.025] \\ \Rightarrow \text{CI OR} &= [e^{0.829},\ e^{1.025}] = [2.29,\ 2.79] \end{align*} \]
Interpretasi Hasil :
Berdasarkan hasil perhitungan, nilai odds ratio bersama (\(\hat{\theta}_{MH}\)) adalah sebesar 2.53. Artinya, setelah mengendalikan variabel \(\textbf{usia}\), individu dengan obesitas memiliki kemungkinan sekitar 2.5 kali lebih besar untuk mengalami diabetes dibandingkan dengan individu yang tidak obesitas.
Interval kepercayaan \(95\%\) antara 2.29 hingga 2.79 tidak mencakup angka 1, sehingga hubungan ini dapat dianggap signifikan secara statistik pada taraf signifikansi \(5\%\).
library(epitools)
odds_ratio_data <- array(
c(1, 9, 2, 8, # Usia <5
5, 9, 6, 8, # Usia 5-14
9, 6, 4, 3, # Usia 15-49
12, 5, 4, 9, # Usia 50-69
14, 4, 5, 9), # Usia >70
dim = c(2, 2, 5),
dimnames = list(
Obesitas = c("Ya", "Tidak"),
Diabetes = c("Ya", "Tidak"),
Usia = c("<5", "5-14", "15-49", "50-69", ">70")
)
)
mh_test <- mantelhaen.test(odds_ratio_data, correct = FALSE)
print(mh_test)
##
## Mantel-Haenszel chi-squared test without continuity correction
##
## data: odds_ratio_data
## Mantel-Haenszel X-squared = 3.8786, df = 1, p-value = 0.04891
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.9919513 4.1659029
## sample estimates:
## common odds ratio
## 2.032824
Dalam analisis data kategorik, khususnya dalam tabel kontingensi tiga arah, sering kali kita ingin mengetahui apakah hubungan antara dua variabel (misalnya X dan Y) konsisten di seluruh strata dari variabel ketiga (misalnya Z). Untuk menjawab pertanyaan ini, kita dapat menggunakan uji homogenitas odds ratio, yang salah satu pendekatannya dilakukan dengan statistik Breslow-Day.
Asosiasi dikatakan homogen jika nilai odds ratio dari setiap strata atau tabel parsial memiliki nilai yang sama. Secara matematis, homogenitas odds ratio antar strata dapat dituliskan sebagai:
\[ \theta_{xy(1)} = \theta_{xy(2)} = \cdots = \theta_{xy(K)} \]
Jika kondisi ini terpenuhi, maka dapat disimpulkan bahwa tidak ada interaksi antara variabel X dan Y terhadap variabel kontrol Z. Artinya, hubungan antara X dan Y tidak bergantung pada level Z.
Namun, jika odds ratio berbeda-beda pada masing-masing strata, maka terdapat indikasi interaksi antara X dan Y yang dipengaruhi oleh Z.
Dalam praktiknya, menguji homogenitas odds ratio sangat penting untuk menentukan apakah kita dapat menyajikan satu nilai odds ratio gabungan (misalnya menggunakan odds ratio dari uji Mantel-Haenszel), atau kita perlu menganalisis tiap strata secara terpisah karena adanya interaksi.
Pengujian homogenitas menggunakan statistik Breslow-Day (BD) dilakukan dengan menetapkan hipotesis sebagai berikut:
Hipotesis Nol (\(H_0\)):
Tidak ada perbedaan odds ratio antar strata, yaitu:
\[ \theta_{xy(1)} = \theta_{xy(2)} = \cdots = \theta_{xy(K)} \]
Hipotesis Alternatif (\(H_1\)):
Setidaknya terdapat satu strata di mana nilai odds ratio berbeda dari yang lain.
Statistik uji yang digunakan untuk menguji homogenitas tersebut adalah sebagai berikut:
\[ X^2_{HBD} = \sum_{j=1}^{K} \frac{(a_j - \tilde{a}_j)^2}{\widehat{\text{Var}}(a_j | \hat{OR}_{MH})} \]
dengan: - \(a_j\) = jumlah kasus terpapar yang diamati pada strata ke-\(j\), - \(\tilde{a}_j\) = jumlah kasus terpapar yang diharapkan berdasarkan odds ratio gabungan \(\hat{OR}_{MH}\), - \(\widehat{\text{Var}}(a_j | \hat{OR}_{MH})\) = varians dari \(\tilde{a}_j\).
\[ df = K - 1 \]
Impletasi ke dalam contoh kasus sebelumnya
\[ OR = \frac{a \times d}{b \times c} \]
Untuk setiap kelompok usia:
\[ \begin{aligned} OR_1 &= \frac{(1 \times 8)}{(9 \times 2)} = \frac{8}{18} = 0.44 \\ OR_2 &= \frac{(5 \times 8)}{(9 \times 6)} = \frac{40}{54} = 0.74 \\ OR_3 &= \frac{(9 \times 6)}{(3 \times 4)} = \frac{54}{12} = 4.5 \\ OR_4 &= \frac{(12 \times 9)}{(5 \times 4)} = \frac{108}{20} = 5.4 \\ OR_5 &= \frac{(14 \times 9)}{(4 \times 5)} = \frac{126}{20} = 6.3 \end{aligned} \]
\[ \hat{OR}_{MH} = \frac{\sum_j \frac{a_j d_j}{N_j}}{\sum_j \frac{b_j c_j}{N_j}} \]
Dengan perhitungan lebih lanjut:
\[ \hat{OR}_{MH} = \frac{0.44 + 0.74 + 4.5 + 5.4 + 6.3}{5} \approx 3.48 \]
\[ -m_{j1} n_{j1} \hat{OR}_{MH} + (n_{j2} - m_{j1} + \hat{OR}_{MH}(m_{j1} + m_{j2})) \bar{a}_j + (1 - \hat{OR}_{MH}) \bar{a}_j^2 = 0 \]
Hasilnya:
\[ \begin{aligned} \bar{a}_1 &\approx 1.5 \\ \bar{a}_2 &\approx 5.6 \\ \bar{a}_3 &\approx 6.8 \\ \bar{a}_4 &\approx 10.4 \\ \bar{a}_5 &\approx 12.3 \end{aligned} \]
\[ \text{Var}(a_j | \hat{OR}_{MH}) = \left( \frac{1}{a_j} + \frac{1}{b_j} + \frac{1}{c_j} + \frac{1}{d_j} \right)^{-1} \]
\[ \begin{aligned} \text{Var}(a_1) &\approx 0.6 \\ \text{Var}(a_2) &\approx 2.5 \\ \text{Var}(a_3) &\approx 3.2 \\ \text{Var}(a_4) &\approx 4.1 \\ \text{Var}(a_5) &\approx 5.3 \end{aligned} \]
\[ X^2_{HBD} = \sum_{j} \frac{(a_j - \bar{a}_j)^2}{\text{Var}(a_j)} \]
\[ X^2_{HBD} = \frac{(1 - 1.5)^2}{0.6} + \frac{(5 - 5.6)^2}{2.5} + \frac{(9 - 6.8)^2}{3.2} + \frac{(12 - 10.4)^2}{4.1} + \frac{(14 - 12.3)^2}{5.3} \approx 4.27 \]
\[ X^2_{HBDT} = X^2_{HBD} \times \frac{\sum_j (a_j - \bar{a}_j)^2}{\sum_j \text{Var}(a_j)} \]
\[ X^2_{HBDT} = 4.27 \times \frac{(0.5^2 + 0.6^2 + 2.2^2 + 1.6^2 + 1.7^2)}{(0.6 + 2.5 + 3.2 + 4.1 + 5.3)} \]
\[ X^2_{HBDT} = 4.27 \times \frac{12.54}{15.7} \approx 3.41 \]
\[ p = 1 - \chi^2(X^2_{HBDT}, df = 1) \]
\[ p = 1 - \chi^2(3.41, df = 1) \approx 0.065 \]
# Membuat data dalam bentuk tabel
data <- data.frame(
Usia = c("<5", "5-14", "15-49", "50-69", ">70"),
a = c(1, 5, 9, 12, 14),
b = c(9, 9, 3, 5, 4),
c = c(2, 6, 4, 4, 5),
d = c(8, 8, 6, 9, 9)
)
# Menghitung OR untuk masing-masing kelompok usia
data$OR <- (data$a * data$d) / (data$b * data$c)
# Menampilkan hasil
print(data[, c("Usia", "OR")])
## Usia OR
## 1 <5 0.4444444
## 2 5-14 0.7407407
## 3 15-49 4.5000000
## 4 50-69 5.4000000
## 5 >70 6.3000000
# Menghitung OR_MH
OR_MH <- sum(data$a * data$d / (data$a + data$b + data$c + data$d)) /
sum(data$b * data$c / (data$a + data$b + data$c + data$d))
# Menampilkan hasil
print(paste("Rasio Odds Mantel-Haenszel: ", round(OR_MH, 2)))
## [1] "Rasio Odds Mantel-Haenszel: 2.53"
# Ekspektasi a_j
data$a_hat <- c(1.5, 5.6, 6.8, 10.4, 12.3)
# Varians a_j
data$Var_a <- c(0.6, 2.5, 3.2, 4.1, 5.3)
# Menampilkan hasil ekspektasi dan varians
print(data[, c("Usia", "a_hat", "Var_a")])
## Usia a_hat Var_a
## 1 <5 1.5 0.6
## 2 5-14 5.6 2.5
## 3 15-49 6.8 3.2
## 4 50-69 10.4 4.1
## 5 >70 12.3 5.3
# Menghitung Statistik X^2_HBD
X2_HBD <- sum((data$a - data$a_hat)^2 / data$Var_a)
# Menampilkan hasil
print(paste("Statistik X^2_HBD: ", round(X2_HBD, 2)))
## [1] "Statistik X^2_HBD: 3.24"
# Koreksi Tarone
tarone_correction <- sum((data$a - data$a_hat)^2) / sum(data$Var_a)
X2_HBDT <- X2_HBD * tarone_correction
# Menampilkan hasil
print(paste("Statistik X^2_HBDT (dengan Koreksi Tarone): ", round(X2_HBDT, 2)))
## [1] "Statistik X^2_HBDT (dengan Koreksi Tarone): 2.25"
# Menghitung p-value
p_value <- 1 - pchisq(X2_HBDT, df = 1) ; p_value
## [1] 0.1334937
Dari hasil perhitungan diperoleh p-value = 0.1333. Berdasarkan aturan keputusan:
Karena p-value = 0.1333 > 0.05, maka gagal menolak \(H_0\). Dengan demikian, odds ratio dianggap seragam di seluruh strata usia, dan tidak ditemukan interaksi yang signifikan antara obesitas dan diabetes setelah dikontrol oleh usia.
Uji Breslow-Day digunakan untuk memeriksa apakah odds ratio seragam di seluruh strata dalam tabel kontingensi tiga arah. Dalam analisis ini, karena odds ratio homogen, maka tidak ada interaksi signifikan antara obesitas dan diabetes setelah mempertimbangkan faktor usia. Oleh karena itu, usia tidak menjadi faktor yang mempengaruhi hubungan antara obesitas dan diabetes dalam dataset ini.
Model Linear Generalisasi (GLM) merupakan perluasan dari model regresi linear klasik yang memungkinkan kita untuk menganalisis data ketika variabel respons tidak mengikuti distribusi normal, atau ketika hubungan antara prediktor dan respons tidak linier secara langsung. Dengan GLM, kita dapat menggunakan pendekatan yang lebih fleksibel melalui pemilihan distribusi dan fungsi link yang sesuai.
GLM terdiri dari tiga komponen utama:
Fungsi Linear Prediktor berupa: \[ \eta = X \beta \]
Suatu distribusi dikatakan termasuk ke dalam exponential family jika fungsi kepadatan probabilitasnya dapat dituliskan dalam bentuk:
\[ f(y; \theta, \phi) = \exp \left\{ \frac{y \theta - b(\theta)}{\phi} + c(y, \phi) \right\} \]
Beberapa distribusi umum yang termasuk dalam exponential family antara lain:
Contoh Distribusi Binomial
Distribusi binomial memiliki fungsi probabilitas sebagai berikut:
\[ P(Y = y) = \binom{n}{y} \pi^y (1 - \pi)^{n - y} \]
Untuk menunjukkan bahwa distribusi ini merupakan bagian dari exponential family, kita ubah bentuknya menjadi:
\[ P(Y = y) = \exp \left\{ \log \binom{n}{y} + y \log \left( \frac{\pi}{1 - \pi} \right) + n \log(1 - \pi) \right\} \]
Dengan substitusi berikut:
Kita dapat menuliskan ulang sebagai bentuk exponential family:
\[ f(y; \theta) = \exp \left\{ y \theta - b(\theta) + \log \binom{n}{y} \right\} \]
Karena sesuai dengan bentuk umum exponential family, maka distribusi binomial memang merupakan bagian dari keluarga ini.
Contoh Distribusi Poisson
Distribusi Poisson memiliki fungsi probabilitas:
\[ P(Y = y) = \frac{\lambda^y e^{-\lambda}}{y!} \]
Kita bisa menuliskannya dalam bentuk exponential family:
\[ P(Y = y) = \exp \left\{ y \log \lambda - \lambda - \log y! \right\} \]
Dengan substitusi:
Maka bentuk exponential family-nya adalah:
\[ f(y; \theta) = \exp \left\{ y \theta - e^{\theta} - \log y! \right\} \]
Distribusi Poisson juga termasuk dalam exponential family.
Regresi logistik merupakan salah satu pendekatan statistik yang digunakan untuk memodelkan hubungan antara satu atau lebih variabel independen dengan variabel dependen biner (diskrit). Berbeda dengan regresi linier yang memberikan keluaran nilai kontinu, regresi logistik menghasilkan keluaran berupa nilai probabilitas antara 0 hingga 1, yang kemudian digunakan untuk klasifikasi ke dalam dua kategori, seperti “ya” atau “tidak”, “sakit” atau “sehat”.
Model ini menggunakan fungsi aktivasi sigmoid untuk mengkonversi output menjadi probabilitas. Fungsi sigmoid memiliki bentuk kurva S (S-shaped), yang dituliskan dalam bentuk rumus:
\[ P(y=1|X) = \frac{1}{1 + e^{-(\beta_0 + \beta_1X_1 + \beta_2X_2 + \dots + \beta_pX_p)}} \]
Di mana:
Regresi logistik merupakan metode yang digunakan untuk menganalisis hubungan antara satu atau lebih variabel independen dengan variabel dependen yang bersifat kategorik. Teknik ini berfungsi untuk mengklasifikasikan data ke dalam kelompok-kelompok diskrit. Umumnya, regresi logistik dimanfaatkan dalam konteks pemodelan prediktif, di mana model bertugas untuk memperkirakan probabilitas suatu entitas termasuk atau tidak termasuk ke dalam suatu kategori tertentu berdasarkan karakteristik inputnya.
Sebagai ilustrasi, nilai 0 sering kali digunakan untuk merepresentasikan kelas negatif, sedangkan nilai 1 mewakili kelas positif. Model ini sangat cocok diterapkan pada kasus klasifikasi biner, yakni ketika variabel target hanya memiliki dua kategori yang saling eksklusif, yaitu 0 dan 1.
Aplikasi Regresi Logistik dalam Dunia Nyata
Regresi logistik digunakan secara luas dalam berbagai bidang, termasuk:
Keunggulan Regresi Logistik
Sangat Cocok untuk Data Biner dan Terpisah secara
Linear
Jika dua kelas dapat dipisahkan secara garis lurus dalam ruang fitur,
regresi logistik bekerja sangat baik. Kelas-kelas yang tidak tumpang
tindih secara signifikan dapat dipisahkan dengan garis batas keputusan
yang linier.
Memberikan Koefisien yang Mudah
Diinterpretasi
Koefisien \(\beta_i\) menunjukkan
perubahan log-odds dari output untuk setiap kenaikan satu unit pada
variabel \(X_i\). Interpretasi
probabilistik ini membuat model lebih intuitif.
Jika kita ingin mengubah log-odds menjadi probabilitas, kita gunakan:
\[ odds = \frac{P}{1-P} \]
dan
\[ \log(odds) = \beta_0 + \beta_1X_1 + \dots + \beta_pX_p \]
Persamaan dan Asumsi dalam Regresi
Logistik
Dalam regresi logistik, prediksi dan probabilitas dihitung menggunakan
sebuah fungsi logistik yang dikenal sebagai fungsi sigmoid. Fungsi ini
memiliki bentuk khas seperti huruf S (S-shaped curve), dan berfungsi
untuk mengubah nilai input yang bersifat kontinu (riil) menjadi nilai
probabilitas dalam rentang antara 0 hingga 1.
Model regresi logistik memanfaatkan nilai keluaran dari fungsi sigmoid sebagai dasar pengambilan keputusan klasifikasi. Apabila nilai probabilitas yang dihasilkan melebihi ambang batas yang telah ditentukan (misalnya 0.5), maka model akan mengklasifikasikan observasi tersebut ke dalam salah satu kategori. Sebaliknya, jika nilai probabilitas lebih kecil dari ambang batas tersebut, maka observasi akan diklasifikasikan ke dalam kategori yang lain.
Fungsi sigmoid merupakan komponen utama dalam regresi logistik yang berfungsi untuk mengkonversi nilai linier menjadi probabilitas antara 0 dan 1. Fungsi sigmoid dirumuskan sebagai:
\[ f(x) = \frac{1}{1 + e^{-x}} \]
Output dari fungsi sigmoid digunakan untuk menentukan klasifikasi. Jika nilai probabilitas melebihi ambang batas (misalnya 0.5), maka prediksi diklasifikasikan sebagai kelas 1 (positif). Jika kurang dari 0.5, maka masuk kelas 0 (negatif).
Contoh Kasus Prediksi Keputusan Kredit
Misalkan kita ingin memprediksi apakah seseorang akan disetujui atau tidak disetujui pengajuan kredit berdasarkan skor keuangan mereka.
Kita akan mensimulasikan data dan menggunakan fungsi sigmoid untuk melihat probabilitas persetujuan.
Simulasi data untuk regresi logistik
# Simulasi data sesuai studi kasus kredit
set.seed(123)
n <- 100
x <- seq(300, 850, length.out = n) # Skor keuangan
# Log odds: intercept -10, slope 0.02
log_odds <- -10 + 0.02 * x
prob <- 1 / (1 + exp(-log_odds)) # Fungsi sigmoid
y <- rbinom(n, 1, prob) # Status disetujui/tidak
# Buat data frame
data <- data.frame(x = x, y = y, prob = prob)
Plot Kurva Sigmoid
# Visualisasi menggunakan base R
plot(x, y, pch = 16, col = "gray60",
xlab = "Skor Keuangan", ylab = "Status / Probabilitas",
main = "Simulasi Persetujuan Kredit dengan Kurva Sigmoid")
lines(x, prob, col = "blue", lwd = 2)
abline(h = 0.5, col = "red", lty = 2)
legend("bottomright",
legend = c("Data Biner (0/1)", "Kurva Sigmoid", "Ambang 0.5"),
col = c("gray60", "blue", "red"),
pch = c(16, NA, NA),
lty = c(NA, 1, 2),
lwd = c(NA, 2, 1),
pt.cex = 1.5,
bty = "n")
Interpretasi:
Bentuk kurva S (sigmoid) menunjukkan bahwa semakin tinggi skor keuangan, semakin besar kemungkinan pengajuan kredit disetujui.
Sekitar skor 500-an adalah titik di mana probabilitas melewati ambang 0.5 (mulai diklasifikasikan sebagai disetujui).
Jika skor keuangan < 500:
Kurva biru berada di bawah 0.5 → model memprediksi probabilitas kecil untuk disetujui.
Sebagian besar titik abu-abu bernilai 0.
Jika skor keuangan > 600:
Kurva sigmoid mendekati nilai 1 → hampir pasti disetujui.
Sebagian besar titik bernilai 1 (disetujui).
Ambang 0.5 sebagai pemisah keputusan klasifikasi efektif dalam menggambarkan transisi dari penolakan ke persetujuan berdasarkan skor.
Kurva sigmoid dalam regresi logistik menunjukkan hubungan non-linear antara variabel prediktor dan probabilitas output. Pendekatan ini efektif untuk klasifikasi biner seperti deteksi penyakit, email spam, dan prediksi ya/tidak.
Fungsi link logit dapat dinyatakan dalam bentuk berikut:
\[ g(\mu) = \log \left( \frac{\mu}{1 - \mu} \right) \]
Model Regresi Logistik
\[ \log \left( \frac{\mu}{1 - \mu} \right) = X\beta \]
Fungsi Inverse Link
\[ \mu = \frac{\exp(X\beta)}{1 + \exp(X\beta)} \]
Estimasi Parameter
Metode estimasi parameter pada model Generalized Linear Model (GLM), termasuk regresi logistik, umumnya menggunakan pendekatan Maximum Likelihood Estimation (MLE).
Fungsi log-likelihood untuk regresi logistik dirumuskan sebagai:
\[ l(\beta) = \sum_{i=1}^{n} \left[ y_i \log(\pi_i) + (1 - y_i) \log(1 - \pi_i) \right] \]
dengan:
\[ \pi_i = \frac{\exp(x_i^T \beta)}{1 + \exp(x_i^T \beta)} \]
Estimasi parameter \(\beta\) dilakukan melalui metode iteratif seperti Newton-Raphson atau algoritma Fisher Scoring.
Contoh Kasus dengan R
# Simulasi data
set.seed(2025)
n <- 200
x <- seq(20, 80, length.out = n) # gunakan x terdistribusi merata
eta <- -6 + 0.15 * x # linear predictor
p <- 1 / (1 + exp(-eta)) # fungsi sigmoid
y <- rbinom(n, 1, prob = p) # output biner
data <- data.frame(x = x, y = y)
# Fit model regresi logistik
model <- glm(y ~ x, data = data, family = binomial)
summary(model)
##
## Call:
## glm(formula = y ~ x, family = binomial, data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -5.95643 0.91146 -6.535 6.36e-11 ***
## x 0.15595 0.02234 6.979 2.97e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 247.64 on 199 degrees of freedom
## Residual deviance: 127.85 on 198 degrees of freedom
## AIC: 131.85
##
## Number of Fisher Scoring iterations: 6
#Visual
# Plot data
plot(data$x, data$y,
pch = 16, col = "gray60",
xlab = "X (Prediktor)", ylab = "Y / Probabilitas",
main = "Kurva Sigmoid Regresi Logistik")
# Prediksi probabilitas model
x_sorted <- sort(data$x)
y_pred <- predict(model, newdata = data.frame(x = x_sorted), type = "response")
# Kurva sigmoid merah
lines(x_sorted, y_pred, col = "red", lwd = 2)
# Garis ambang 0.5
abline(h = 0.5, col = "blue", lty = 2)
# Legend lengkap
legend("topleft",
legend = c("Data Biner (0/1)", "Kurva Sigmoid", "Ambang 0.5"),
col = c("gray60", "red", "blue"),
pch = c(16, NA, NA),
lty = c(NA, 1, 2),
lwd = c(NA, 2, 1),
pt.cex = 1.5,
bty = "n")
Interpretasi Output
Koefesien Model
Intercept (β₀): -5.95643
Artinya, ketika nilai prediktor x
= 0, log-odds dari
keberhasilan (y = 1
) adalah -5.95643. Namun, karena
x
pada simulasi ini dimulai dari 20, interpretasi pada
x = 0
kurang relevan secara praktikal.
x (β₁): 0.15595
Setiap kenaikan satu unit pada x
akan meningkatkan log-odds
keberhasilan sebesar 0.15595. Dalam bentuk odds rasio:
exp(0.15595) ≈ 1.169
Artinya, setiap kenaikan satu unit pada x
meningkatkan
peluang y = 1
sekitar 16.9%
Statistik Uji
Nilai p untuk kedua koefisien sangat kecil
(< 0.001
), yang berarti kedua variabel sangat signifikan
dalam model ini.
AIC = 131.85, digunakan untuk membandingkan model; semakin kecil semakin baik.
Iterasi algoritma Fisher Scoring selesai dalam 6 langkah, menunjukkan konvergensi yang cepat.
Model regresi Poisson digunakan saat variabel respons berbentuk data cacah (count data), yaitu bilangan bulat non-negatif seperti jumlah kasus, frekuensi kejadian, dan sebagainya. Model ini merupakan bagian dari Generalized Linear Model (GLM), di mana asumsi utamanya adalah distribusi dari variabel respons mengikuti distribusi Poisson.
Distribusi Poisson memiliki bentuk fungsi probabilitas sebagai berikut:
\[ P(Y = y) = \frac{e^{-\lambda} \lambda^y}{y!} \]
dengan:
Distribusi Poisson termasuk dalam keluarga eksponensial (exponential family) yang umum digunakan dalam GLM.
Distribusi Poisson dapat ditulis dalam bentuk exponential family sebagai:
\[ f(y; \theta) = \exp \left\{ y \log(\lambda) - \lambda - \log(y!) \right\} \]
dengan parameterisasi:
Karena distribusi ini bisa dinyatakan dalam bentuk tersebut, maka ia memenuhi syarat sebagai anggota keluarga eksponensial.
Fungsi Link Dalam konteks GLM, fungsi link digunakan untuk menghubungkan nilai harapan \(\mu = E(Y)\) dengan prediktor linear.
Fungsi link kanonik untuk distribusi Poisson adalah:
\[ g(\mu) = \log(\mu) \]
Maka model regresinya menjadi:
\[ \log(\mu_i) = \mathbf{x}_i^\top \boldsymbol{\beta} \]
dengan:
Fungsi inversenya:
\[ \mu_i = \exp(\mathbf{x}_i^\top \boldsymbol{\beta}) \]
Fungsi eksponensial menjamin bahwa nilai \(\mu_i\) selalu positif.
Estimasi Parameter
Parameter \(\boldsymbol{\beta}\) diestimasi menggunakan metode Maximum Likelihood Estimation (MLE). Fungsi log-likelihood dari model regresi Poisson:
\[ \ell(\boldsymbol{\beta}) = \sum_{i=1}^{n} \left[ y_i (\mathbf{x}_i^\top \boldsymbol{\beta}) - \exp(\mathbf{x}_i^\top \boldsymbol{\beta}) - \log(y_i!) \right] \]
Estimasi dapat dilakukan melalui metode numerik seperti:
Kelebihan Model Regresi Poisson - Cocok untuk data cacah - Termasuk dalam kerangka GLM - Parameter mudah diinterpretasi dalam skala log atau eksponensial
Keterbatasan Model Regresi Poisson
Jika terjadi overdispersion, bisa digunakan model alternatif seperti:
Contoh Kasus misal kita memiliki data simulasi sebagai berikut
#data
set.seed(123)
n <- 150
x <- runif(n, 0, 10) # variabel prediktor (misal: waktu aktivitas fisik per minggu)
lambda <- exp(1 + 0.3 * x) # hubungan log-linear
y <- rpois(n, lambda) # respon (jumlah keluhan kesehatan per bulan)
data <- data.frame(y, x)
Estimasi Regresi Poisson
poisson_model <- glm(y ~ x, data = data, family = poisson)
summary(poisson_model)
##
## Call:
## glm(formula = y ~ x, family = poisson, data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.014638 0.062238 16.3 <2e-16 ***
## x 0.298834 0.008187 36.5 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 1738.8 on 149 degrees of freedom
## Residual deviance: 127.4 on 148 degrees of freedom
## AIC: 778.6
##
## Number of Fisher Scoring iterations: 4
Plot Hasil Prediksi
plot(data$x, data$y, pch = 16, col = "darkgray", main = "Data dan Hasil Prediksi")
newdata <- data.frame(x = seq(min(x), max(x), length.out = 100))
pred <- predict(poisson_model, newdata = newdata, type = "response")
lines(newdata$x, pred, col = "blue", lwd = 2)
Diagnostik dan Overdispersion
dispersion <- sum(residuals(poisson_model, type = "pearson")^2) / poisson_model$df.residual
dispersion
## [1] 0.8601197
Interpretasi : Model regresi Poisson yang diestimasi memiliki bentuk:
\[ \log(\mu_i) = 1.014638 + 0.298834 \cdot x_i \]
atau jika dinyatakan dalam bentuk eksponensial:
\[ \mu_i = \exp(1.014638 + 0.298834 \cdot x_i) \]
Koefisien intercept sebesar 1.014638 menyatakan bahwa ketika \(x = 0\), maka rata-rata jumlah kejadian (\(\mu\)) yang diprediksi adalah sebesar \(\exp(1.014638) \approx 2.76\). Sementara itu, koefisien untuk variabel \(x\) sebesar 0.298834 menunjukkan bahwa setiap peningkatan satu satuan pada \(x\), akan meningkatkan rata-rata jumlah kejadian (\(\mu\)) sebesar \(\exp(0.298834) \approx 1.35\) kali, atau sekitar 35% peningkatan.
Secara statistik, kedua koefisien memiliki nilai p yang sangat kecil (\(< 2 \times 10^{-16}\)), menunjukkan bahwa keduanya signifikan dalam memprediksi variabel respons. Penurunan deviance dari nilai null deviance sebesar 1738.8 menjadi residual deviance sebesar 127.4 menandakan bahwa penambahan variabel prediktor \(x\) secara substansial meningkatkan kecocokan model. Nilai AIC yang relatif rendah (778.6) juga mengindikasikan bahwa model memiliki performa yang baik dalam menjelaskan data.
Untuk menguji keberadaan overdispersion, dilakukan perhitungan menggunakan residual Pearson dengan rumus:
\[ \text{Dispersion} = \frac{\sum (\text{residual Pearson})^2}{\text{df residual}} \]
Hasilnya diperoleh nilai dispersion sebesar 0.8601197. Karena nilai ini kurang dari 1, maka dapat disimpulkan bahwa tidak terdapat overdispersion dalam model. Dengan kata lain, asumsi distribusi Poisson yang mengharuskan bahwa varians sama dengan nilai harapan (mean) telah terpenuhi. Oleh karena itu, model Poisson ini dapat dianggap layak untuk digunakan dalam analisis lebih lanjut, baik untuk keperluan inferensi maupun prediksi.
Dalam analisis Generalized Linear Model (GLM), salah satu aspek penting dalam inferensi statistik adalah memahami ekspektasi dan varians dari estimator parameter. Pemahaman ini sangat diperlukan dalam pengembangan alat uji statistik seperti uji Wald, uji Likelihood Ratio, serta dalam membentuk interval kepercayaan.
Ekspektasi dan Varians dalam GLM
\[ \mathbb{E}[\hat{\beta}] = \beta \]
Dalam konteks GLM, estimasi parameter biasanya diperoleh dengan metode Maximum Likelihood Estimation (MLE). Estimator MLE bersifat asymptotically unbiased, yang berarti bahwa seiring bertambahnya jumlah sampel, bias dari estimator semakin kecil hingga mendekati nol.
Varians mengukur presisi dari estimator. Semakin kecil varians, semakin presisi estimasi tersebut. Secara umum, varians dari estimator \(\hat{\beta}\) dalam GLM didekati oleh:
\[ \text{Var}(\hat{\beta}) \approx (\mathbf{X}^\top \mathbf{W} \mathbf{X})^{-1} \]
Di mana:
Distribusi Asimptotik Estimator Untuk ukuran sampel yang besar, distribusi dari estimator \(\hat{\beta}\) mendekati distribusi normal:
\[ \hat{\beta} \sim \mathcal{N}(\beta, \text{Var}(\hat{\beta})) \]
Distribusi ini digunakan sebagai dasar dari:
Varians dalam GLM Tidak Konstan
Berbeda dengan regresi linear klasik (OLS) yang mengasumsikan homoskedastisitas dengan bentuk:
\[ \text{Var}(Y_i) = \sigma^2 \]
GLM mengizinkan heteroskedastisitas, dan varians data bergantung pada nilai harapan (mean) dari respon. Secara umum dalam GLM:
\[ \text{Var}(Y_i) = \phi V(\mu_i) \]
Dimana:
Beberapa contoh fungsi varians dalam GLM:
Mengapa Asumsi Ini Penting?
Dalam GLM, keakuratan dari standar error, p-value, dan interval kepercayaan sangat bergantung pada validitas asumsi mengenai bentuk varians dan distribusi. Jika model salah spesifikasi (misalnya menggunakan distribusi Poisson untuk data yang overdispersed), maka inferensi statistik bisa menjadi menyesatkan. Oleh karena itu, diagnostik model dan pemeriksaan asumsi sangat penting dalam praktik analisis GLM.
Selain itu, keberadaan parameter dispersi \(\phi\) membuka peluang untuk mengukur overdispersi dalam data. Jika \(\phi > 1\), data dianggap memiliki overdispersi, yang dapat diatasi dengan menggunakan model alternatif seperti Quasi-Poisson atau Negative Binomial.
Contoh Kasus Misalkan sebuah dinas transportasi ingin menganalisis hubungan antara intensitas lalu lintas (jumlah kendaraan per jam) dan jumlah kecelakaan yang terjadi di suatu ruas jalan dalam sehari. Karena data berupa hitungan (jumlah kecelakaan), maka digunakan model regresi Poisson.
Simulasi data dan Estimasi Model
# Simulasi data (dibuat beda dari contoh lain!)
set.seed(321)
traffic <- runif(120, min = 50, max = 500) # jumlah kendaraan per jam
lambda <- exp(1 + 0.004 * traffic) # rata-rata kecelakaan
accidents <- rpois(120, lambda)
# Estimasi model regresi Poisson
model2 <- glm(accidents ~ traffic, family = poisson)
# Ringkasan hasil model
summary(model2)
##
## Call:
## glm(formula = accidents ~ traffic, family = poisson)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.0289841 0.0892996 11.52 <2e-16 ***
## traffic 0.0038906 0.0002453 15.86 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 378.34 on 119 degrees of freedom
## Residual deviance: 107.30 on 118 degrees of freedom
## AIC: 579.69
##
## Number of Fisher Scoring iterations: 4
Visualisasi plot
# Visualisasi
library(ggplot2)
ggplot(data = data.frame(traffic, accidents), aes(x = traffic, y = accidents)) +
geom_point(color = "darkred", alpha = 0.7) +
stat_smooth(method = "glm", method.args = list(family = "poisson"),
se = FALSE, color = "navyblue", size = 1) +
labs(
title = "Hubungan Antara Lalu Lintas dan Jumlah Kecelakaan",
x = "Jumlah Kendaraan per Jam",
y = "Jumlah Kecelakaan per Hari"
) +
theme_minimal(base_size = 13) +
theme(
plot.title = element_text(face = "bold", hjust = 0.5),
panel.grid.major = element_line(color = "gray80"),
panel.grid.minor = element_blank()
)
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## `geom_smooth()` using formula = 'y ~ x'
Berdasarkan hasil estimasi model regresi Poisson dengan variabel
respon accidents
dan prediktor traffic
,
diperoleh model sebagai berikut:
\[ \log(\hat{\mu}_i) = \hat{\beta}_0 + \hat{\beta}_1 \cdot \text{traffic}_i \]
dengan nilai estimasi parameter sebesar: - Intersep \(\hat{\beta}_0 = 1.029\) - Koefisien trafik \(\hat{\beta}_1 = 0.00389\)
Kedua parameter tersebut signifikan secara statistik dengan nilai p-value < 0.001, yang menunjukkan bahwa hubungan antara jumlah kendaraan dan jumlah kecelakaan sangat signifikan secara statistik.
Koefisien trafik sebesar 0.00389 dapat diinterpretasikan sebagai peningkatan log-rata-rata jumlah kecelakaan untuk setiap kenaikan 1 unit kendaraan per jam. Jika dikonversi ke skala asli melalui eksponensial, diperoleh:
\[ \exp(0.00389) \approx 1.0039 \]
yang berarti, setiap kenaikan 1 kendaraan per jam diperkirakan akan meningkatkan jumlah kecelakaan sebesar 0.39%.
Model memiliki null deviance sebesar 378.34 dan residual deviance sebesar 107.30 dengan derajat bebas masing-masing 119 dan 118. Hal ini menunjukkan bahwa model berhasil menjelaskan sebagian besar variasi dalam data. Nilai Akaike Information Criterion (AIC) sebesar 579.69 menunjukkan model yang cukup baik, mengingat semakin kecil nilai AIC, semakin baik kualitas model.
Asumsi distribusi Poisson masih terpenuhi karena nilai dispersion parameter adalah 1. Selain itu, rasio residual deviance terhadap derajat bebas adalah:
\[ \frac{107.30}{118} \approx 0.91 \]
Nilai ini mendekati 1, sehingga tidak ditemukan indikasi kuat adanya overdispersi. Dengan demikian, tidak diperlukan model alternatif seperti quasipoisson atau negative binomial.
Penyesuaian Konsep Ekspektasi dan Varians
Sesuai dengan konsep Generalized Linear Model (GLM), pada model Poisson:
Parameter \(\hat{\beta}\) diperoleh menggunakan metode Maximum Likelihood Estimation (MLE), yang dalam GLM bersifat asimtotik tak bias (asymptotically unbiased) dan memiliki varians:
\[ \text{Var}(\hat{\beta}) \approx (\mathbf{X}^\top \mathbf{W} \mathbf{X})^{-1} \]
Distribusi asimptotik dari estimator mengikuti:
\[ \hat{\beta} \sim \mathcal{N}(\beta, \text{Var}(\hat{\beta})) \]
Distribusi ini merupakan dasar dari pengujian statistik seperti uji Wald, pembentukan interval kepercayaan, dan penghitungan nilai p-value.
Kesimpulan
Ekspektasi Dalam kerangka Generalized Linear Model (GLM), ekspektasi dari variabel respon dapat diturunkan secara teoritis dari bentuk distribusi eksponensial. Distribusi eksponensial umum dituliskan dalam bentuk fungsi kepadatan sebagai:
\[ f(y; \theta) = \exp\{y \theta - b(\theta) + c(y)\} \]
Fungsi log-likelihood dari distribusi eksponensial ini adalah:
\[ \log f(y; \theta) = y \theta - b(\theta) + c(y) \]
Berdasarkan teori fungsi momen, ekspektasi dari \(Y\) diperoleh dari:
\[ \mathbb{E}(Y) = \int y f(y;\theta) dy = \mu \]
Untuk mendapatkan nilai ekspektasi secara matematis, kita menggunakan turunan pertama dari log-likelihood terhadap parameter \(\theta\), yang dikenal sebagai fungsi skor:
\[ U(\theta) = \frac{\partial \ell}{\partial \theta} = y - b'(\theta) \]
Mengambil ekspektasi dari fungsi skor, dan mengingat bahwa MLE adalah solusi dari nilai ekspektasi turunan pertama sama dengan nol:
\[ \mathbb{E}[U(\theta)] = \mathbb{E}[y - b'(\theta)] = \mu - b'(\theta) = 0 \]
Maka dapat disimpulkan bahwa:
\[ \mu = b'(\theta) \]
Artinya, ekspektasi dari variabel respon dalam GLM bergantung pada bentuk turunan pertama dari fungsi \(b(\theta)\), yaitu fungsi kumulatif dari parameter kanonik \(\theta\).
Varians Selanjutnya, untuk menurunkan varians dari variabel respon \(Y\), kita menggunakan turunan kedua dari fungsi log-likelihood:
\[ \frac{\partial^2 \ell}{\partial \theta^2} = -b''(\theta) \]
Varians dari estimator MLE diberikan oleh bentuk informasi Fisher:
\[ \text{Var}(Y) = \left[ -\mathbb{E} \left( \frac{\partial^2 \ell}{\partial \theta^2} \right) \right]^{-1} \]
Namun dalam GLM, secara langsung diperoleh bahwa:
\[ \text{Var}(Y) = b''(\theta) \]
Dalam bentuk umum GLM yang memasukkan faktor dispersi \(\phi\), maka varian diperoleh sebagai:
\[ \text{Var}(Y) = \phi V(\mu) \]
di mana \(V(\mu) = b''(\theta)\) adalah fungsi varian yang tergantung pada distribusi dari data.
Fungsi Varians Berdasarkan Jenis Distribusi Beberapa distribusi umum yang termasuk dalam keluarga eksponensial memiliki bentuk fungsi varian \(V(\mu)\) yang khas:
Pemahaman mengenai bentuk fungsi varian ini penting karena akan mempengaruhi bentuk matriks kovarians dari parameter \(\hat{\beta}\), serta memengaruhi ketepatan dan presisi inferensi yang dilakukan menggunakan GLM.
Maximum Likelihood Estimation (MLE) Maximum Likelihood Estimation (MLE) merupakan pendekatan standar dalam statistika untuk memperkirakan parameter model. Tujuan utamanya adalah mencari nilai parameter yang memaksimalkan fungsi likelihood, yaitu probabilitas dari data yang diamati diberikan parameter tertentu.
Konsep dasar: Pendekatan ini didasarkan pada ide bahwa nilai parameter terbaik adalah yang memaksimalkan nilai dari fungsi likelihood, atau secara lebih umum, log-likelihood.
Langkah-langkah utama dalam MLE:
Menentukan titik stasioner dari fungsi log-likelihood, yaitu dengan menset: \[ \frac{\partial \ell(\beta)}{\partial \beta} = 0 \]
Memastikan bahwa titik tersebut adalah maksimum dengan melihat tanda dari turunan keduanya: \[ \frac{\partial^2 \ell(\beta)}{\partial \beta^2} < 0 \]
Tantangan dalam GLM Dalam konteks Generalized Linear Models (GLM), bentuk eksplisit dari fungsi log-likelihood tidak selalu tersedia, atau terlalu kompleks untuk diselesaikan secara analitik. Oleh karena itu, digunakan metode numerik untuk mendekati solusi dari estimasi parameter.
Metode Optimisasi Newton-Raphson
Newton-Raphson adalah metode iteratif yang digunakan untuk menyelesaikan sistem nonlinear, termasuk dalam konteks MLE.
Karakteristik utama: - Menggunakan vektor gradien (score vector) dari fungsi log-likelihood. - Menggunakan matriks Hessian, yaitu matriks turunan kedua dari fungsi log-likelihood.
Formula Iteratif Newton-Raphson:
\[ \beta^{(t+1)} = \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \]
Di mana: - \(\beta^{(t)}\): nilai parameter pada iterasi ke-\(t\) - \(U(\beta^{(t)})\): vektor score, yaitu turunan pertama log-likelihood - \(H(\beta^{(t)})\): matriks Hessian, yaitu turunan kedua log-likelihood
Fisher Scoring Fisher Scoring adalah modifikasi dari Newton-Raphson, di mana matriks Hessian \(H(\beta)\) digantikan dengan matriks informasi Fisher \(I(\beta)\). Hal ini dilakukan karena matriks informasi Fisher memiliki ekspektasi yang lebih stabil.
\[ \beta^{(t+1)} = \beta^{(t)} + I^{-1}(\beta^{(t)}) U(\beta^{(t)}) \]
IRLS (Iteratively Reweighted Least Squares)
Merupakan pendekatan numerik lain yang juga merupakan modifikasi dari Fisher Scoring. IRLS digunakan dalam GLM karena memungkinkan penggunaan metode kuadrat terkecil dengan pembobotan ulang di setiap iterasi, sehingga estimasi parameter menjadi mirip dengan regresi biasa, tetapi lebih fleksibel untuk berbagai distribusi data.
Implementasi Newton-Raphson
Komponen-komponen penting:
Statistik Score ke-\(j\): \[ U_j(\beta) = \frac{\partial \log L(\beta)}{\partial \beta_j} \]
Matriks Hessian: \[ H_{jk}(\beta) = \frac{\partial^2 \log L(\beta)}{\partial \beta_j \partial \beta_k} \] Ekspansi Taylor:
Dalam pendekatan numerik, digunakan pendekatan deret Taylor terhadap fungsi score:
\[ U(\beta^*) \approx U(\beta) + H(\beta)(\beta^* - \beta) \]
Estimasi Parameter (Iteratif):
Dari ekspansi di atas, diperoleh estimasi parameter secara iteratif:
\[ \hat{\beta} \approx \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \]
Catatan tambahan: - Metode Newton-Raphson dan turunannya sangat penting dalam GLM karena distribusi dari data seringkali kompleks dan tidak memungkinkan solusi eksplisit. - Kombinasi antara teori turunan (differential calculus), probabilitas, dan optimasi menjadi pondasi utama dari proses estimasi dalam GLM.
Diagnostik dalam Generalized Linear Model (GLM) bertujuan untuk mengevaluasi sejauh mana model yang telah dibangun mampu merepresentasikan data dengan baik. Proses ini mencakup pemeriksaan kesesuaian model (goodness of fit), deteksi adanya outlier atau data pencilan, serta evaluasi terhadap asumsi-asumsi model.
Beberapa pendekatan yang umum digunakan untuk evaluasi model GLM antara lain:
Statistik Devians
Statistik devians merupakan ukuran kesesuaian model yang membandingkan model terestimasi terhadap model saturated (model yang sangat kompleks dan memiliki kesesuaian sempurna terhadap data).
Rumus devians dinyatakan sebagai:
\[ D = 2 \sum \left[ y_i \log \left( \frac{y_i}{\hat{\mu}_i} \right) - (y_i - \hat{\mu}_i) \right] \]
Keterangan: - \(y_i\) adalah nilai aktual, - \(\hat{\mu}_i\) adalah nilai prediksi dari model, - Model saturated memberikan devians minimum (nol), sedangkan devians dari model GLM dibandingkan terhadapnya.
Statistik Chi-Kuadrat Pearson Statistik ini digunakan untuk menguji apakah model memberikan hasil yang lebih baik dibandingkan dengan model tanpa penjelas sama sekali (null model).
\[ X^2 = \sum \frac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i} \]
Interpretasi: - Nilai \(X^2\) yang besar mengindikasikan ketidaksesuaian model. - Jika hasil uji signifikan (misalnya p-value < 0.05), maka dapat disimpulkan bahwa model lebih baik daripada tidak ada model sama sekali.
Catatan Tambahan: - Jika data dikelompokkan, maka statistik devians dan statistik Chi-Kuadrat Pearson mengikuti distribusi Chi-Square. Ini memungkinkan dilakukannya uji hipotesis formal. - Jika data tidak dikelompokkan (misalnya data individual), distribusi Chi-Square tidak lagi berlaku, dan pendekatan lain diperlukan untuk evaluasi. - Karena MLE (Maximum Likelihood Estimation) meminimalkan devians, statistik devians menjadi sangat berguna untuk evaluasi model berbasis MLE.
Analisis Residual
Residual merupakan selisih antara nilai observasi dan nilai yang diprediksi oleh model:
\[ e_i = y_i - \hat{\mu}_i \]
Manfaat analisis residual: - Mengidentifikasi pola yang tidak ditangkap oleh model. - Mendeteksi penyimpangan sistematis, seperti heteroskedastisitas atau non-linearitas. - Membantu mengevaluasi asumsi-asumsi model, terutama melalui plot residual vs fitted value.
Plot residual yang ideal tidak menunjukkan pola tertentu (acak) — jika terdapat pola sistematis, maka bisa jadi model belum tepat atau ada asumsi yang dilanggar.
Regresi logistik merupakan salah satu jenis model regresi yang digunakan untuk memodelkan hubungan antara satu atau lebih variabel prediktor (bebas) dengan variabel respons yang bersifat biner (memiliki dua kategori, misalnya 0 dan 1). Model ini sangat populer dalam analisis data kategorik, terutama ketika ingin memprediksi probabilitas terjadinya suatu peristiwa.
Tujuan utama dari regresi logistik adalah untuk memperkirakan probabilitas bahwa suatu observasi masuk ke dalam salah satu dari dua kategori. Oleh karena itu, bentuk modelnya bersifat non-linear terhadap parameter, dan proses estimasi dilakukan melalui pendekatan Maximum Likelihood Estimation (MLE).
Model regresi logistik secara umum dinyatakan dalam bentuk:
\[ \pi(x) = \frac{\exp(\beta_0 + \beta_1 x)}{1 + \exp(\beta_0 + \beta_1 x)} \]
Fungsi ini disebut sebagai fungsi logistik, yang membatasi hasil prediksi berada dalam rentang antara 0 hingga 1.
Untuk \(n\) observasi, fungsi log-likelihood dari model regresi logistik dinyatakan sebagai:
\[ \ell(\beta) = \sum_{i=1}^{n} \left[y_i \log(\pi_i) + (1 - y_i)\log(1 - \pi_i)\right] \]
di mana: - \(y_i\): nilai aktual dari respons, - \(\pi_i\): probabilitas prediksi dari model untuk observasi ke-\(i\).
Estimasi dengan Newton-Raphson Metode Newton-Raphson digunakan untuk menemukan nilai parameter \(\beta\) yang memaksimalkan fungsi log-likelihood tersebut. Metode ini termasuk dalam kelas metode numerik, yang melibatkan kalkulasi turunan pertama dan kedua.
Untuk observasi ke-\(i\), fungsi probabilitasnya adalah:
\[ \pi_i = \frac{1}{1 + \exp(-\mathbf{x}_i^\top \beta)} \]
Sehingga log-likelihood keseluruhan menjadi:
\[ \ell(\beta) = \sum_{i=1}^{n} \left[y_i \log(\pi_i) + (1 - y_i)\log(1 - \pi_i)\right] \] Langkah-Langkah Newton-Raphson
Prosedur Newton-Raphson untuk regresi logistik melibatkan tiga komponen utama:
Score function merupakan turunan pertama dari fungsi log-likelihood terhadap parameter:
\[ \mathbf{U}(\beta) = \frac{\partial \ell(\beta)}{\partial \beta} = \mathbf{X}^\top (\mathbf{y} - \boldsymbol{\pi}) \]
di mana: - \(\mathbf{X}\): matriks desain (n x p), - \(\mathbf{y}\): vektor nilai observasi, - \(\boldsymbol{\pi}\): vektor probabilitas prediksi.
Matriks Hessian adalah turunan kedua dari log-likelihood:
\[ \mathbf{H}(\beta) = -\mathbf{X}^\top \mathbf{W} \mathbf{X} \]
dengan: - \(\mathbf{W} = \text{diag}(\pi_i(1 - \pi_i))\), sebuah matriks diagonal yang berisi variansi dari distribusi Bernoulli.
Langkah iteratif Newton-Raphson untuk memperbaharui estimasi parameter dinyatakan sebagai:
\[ \beta^{(t+1)} = \beta^{(t)} + \left(\mathbf{X}^\top \mathbf{W}^{(t)} \mathbf{X} \right)^{-1} \mathbf{X}^\top (\mathbf{y} - \boldsymbol{\pi}^{(t)}) \]
Proses ini dilakukan secara berulang hingga konvergen, yaitu ketika perubahan nilai \(\beta\) antar iterasi sudah sangat kecil atau berada di bawah ambang batas tertentu.
Estimasi MLE dengan Newton-Raphson (Manual di R)
# Set seed untuk replikasi
set.seed(123)
# Jumlah data
n <- 150
# Variabel prediktor (dua fitur)
x1 <- rnorm(n, mean = 2, sd = 1.5)
x2 <- runif(n, min = -1, max = 3)
# Gabungkan ke dalam matriks X dengan intercept
X <- cbind(1, x1, x2)
# Parameter sebenarnya (true beta)
beta_true <- c(-0.8, 1.2, -1.5)
# Hitung eta dan probabilitas
eta <- X %*% beta_true
p <- 1 / (1 + exp(-eta))
# Buat respon biner dari distribusi Bernoulli
y <- rbinom(n, size = 1, prob = p)
Newton-Raphson Iterasi Manual
# Inisialisasi
beta <- rep(0, ncol(X))
tol <- 1e-6
max_iter <- 100
for (i in 1:max_iter) {
eta <- X %*% beta
pi_hat <- 1 / (1 + exp(-eta))
W <- diag(as.numeric(pi_hat * (1 - pi_hat)))
z <- eta + solve(W) %*% (y - pi_hat)
beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z
if (sum(abs(beta_new - beta)) < tol) break
beta <- beta_new
}
# Tampilkan hasil estimasi
beta
## [,1]
## -1.461466
## x1 1.410017
## x2 -1.389667
Kesimpulan Estimasi Estimasi parameter 𝛽 dalam model regresi logistik dilakukan melalui metode Maximum Likelihood Estimation (MLE).
Proses estimasi ini menggunakan algoritma Newton-Raphson, yaitu metode numerik iteratif untuk menemukan nilai maksimum dari fungsi log-likelihood.
Newton-Raphson menggunakan pendekatan matematis berupa:
Turunan pertama (score function) untuk menentukan arah perbaikan,
Turunan kedua (Hessian matrix) untuk menyesuaikan kecepatan konvergensi.
Dalam kasus ini, kita menggunakan dua variabel prediktor yang disimulasikan secara acak dari distribusi normal dan uniform untuk memastikan struktur data berbeda dari contoh standar.
Hasil estimasi mendekati nilai parameter sebenarnya yang digunakan dalam simulasi, menunjukkan bahwa algoritma bekerja dengan baik.
Inferensi Parameter
Uji Wald digunakan untuk mengevaluasi apakah suatu parameter regresi logistik (\(\beta_j\)) signifikan secara statistik atau tidak. Artinya, kita ingin mengetahui apakah variabel prediktor ke-\(j\) memiliki kontribusi yang nyata terhadap model.
Hipotesis yang diuji adalah:
\[ \begin{aligned} H_0 &: \beta_j = 0 \quad \text{(tidak signifikan)} \\ H_1 &: \beta_j \ne 0 \quad \text{(signifikan)} \end{aligned} \]
Jika hipotesis nol diterima, maka dapat dikatakan bahwa tidak terdapat bukti yang cukup untuk menyatakan bahwa \(\beta_j\) berbeda dari nol.
Teori Dasar Uji Wald
Uji ini didasarkan pada teori bahwa estimator Maximum Likelihood Estimation (MLE) untuk parameter \(\beta_j\) memiliki distribusi mendekati normal ketika ukuran sampel cukup besar, yaitu:
\[ \hat{\beta}_j \sim \mathcal{N}(\beta_j, \text{Var}(\hat{\beta}_j)) \]
Apabila \(H_0: \beta_j = 0\) benar, maka statistik uji Wald dapat dihitung dengan menggunakan:
\[ Z = \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)} \sim \mathcal{N}(0, 1) \]
Untuk mempermudah interpretasi, Z dikuadratkan dan diperoleh statistik Wald sebagai berikut:
\[ W = Z^2 = \left( \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)} \right)^2 \sim \chi^2_1 \]
Nilai statistik Wald ini kemudian dibandingkan dengan nilai kritis dari distribusi chi-square derajat kebebasan 1, atau dapat juga digunakan untuk menghitung nilai-p (p-value) dalam pengambilan keputusan.
Interpretasi Uji Wald
Simulasi dan Uji Wald Langkah demi Langkah:
# Simulasi data yang berbeda
set.seed(42)
n <- 120
x1 <- rnorm(n, mean = 1)
x2 <- runif(n, -1, 2)
X <- cbind(1, x1, x2) # tambahkan intercept
beta_true <- c(-1.2, 0.7, 1.5)
eta <- X %*% beta_true
p <- 1 / (1 + exp(-eta))
y <- rbinom(n, 1, p)
# Masukkan data simulasi ke dalam data frame
data <- data.frame(y = y, x1 = x1, x2 = x2)
# Fit model regresi logistik
model <- glm(y ~ x1 + x2, data = data, family = binomial)
summary(model)
##
## Call:
## glm(formula = y ~ x1 + x2, family = binomial, data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.0786 0.3538 -3.049 0.00230 **
## x1 0.7836 0.2402 3.262 0.00111 **
## x2 1.4994 0.3111 4.820 1.43e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 165.52 on 119 degrees of freedom
## Residual deviance: 122.17 on 117 degrees of freedom
## AIC: 128.17
##
## Number of Fisher Scoring iterations: 4
Langkah 1: Ambil Nilai Koefesien dan Standar Error
beta_hat <- coef(model)["x2"]
se_beta <- summary(model)$coefficients["x2", "Std. Error"]
Langkah 2: Hitung Statistik Z
Z <- beta_hat / se_beta
Z
## x2
## 4.820372
Langkah 3: Hitung Statistik Wald
Wald_stat <- Z^2
Wald_stat
## x2
## 23.23599
Langkah 4: Hitung p-value
p_value <- 1 - pchisq(Wald_stat, df = 1)
p_value
## x2
## 1.432904e-06
Jika p-value < 0.05, maka parameter 𝛽 β signifikan, artinya variabel prediktor secara statistik berpengaruh terhadap peluang keberhasilan.
Jika p-value > 0.05, maka tidak terdapat cukup bukti untuk menyatakan bahwa variabel tersebut berpengaruh.
Kesimpulan
Berdasarkan hasil perhitungan Uji Wald pada model regresi logistik, diperoleh nilai statistik Z sebesar 4.820 dan statistik Wald sebesar 23.236 dengan p-value yang sangat kecil yaitu 1.43×10⁻⁶. Karena p-value < 0.05, maka dapat disimpulkan bahwa koefisien parameter untuk variabel x2 signifikan secara statistik. Artinya, terdapat cukup bukti untuk menolak hipotesis nol, sehingga variabel x2 berpengaruh signifikan terhadap peluang terjadinya kejadian dalam model regresi logistik yang digunakan
Banding model penuh dengan model tanpa prediktor
# Model null
model_null <- glm(y ~ 1, data = data, family = binomial)
# Likelihood ratio test
anova(model_null, model, test = "Chisq")
## Analysis of Deviance Table
##
## Model 1: y ~ 1
## Model 2: y ~ x1 + x2
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 119 165.52
## 2 117 122.17 2 43.347 3.866e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Evaluasi Kebaikan Model
Akaike Information Criterion (AIC)
Semakin Kecil AIC, semakin baik model
AIC(model)
## [1] 128.174
Bayesian Information Criterion(BIC)
Alternatif terhadap AIC, menghukum kompleksitas model.
BIC(model)
## [1] 136.5365
Kesimpulan
Estimasi regresi logistik dilakukan dengan metode Maximum Likelihood Estimation (MLE) menggunakan iterasi Newton-Raphson.
Berdasarkan uji likelihood ratio, nilai deviance sebesar 43.347 dengan p-value sebesar 3.866e-10, yang jauh lebih kecil dari 0.05. Hal ini menunjukkan bahwa model dengan prediktor (x1 dan x2) secara signifikan lebih baik daripada model null (tanpa prediktor).
Nilai AIC sebesar 128.174 dan BIC sebesar 136.5365, menunjukkan model memiliki keseimbangan yang baik antara kecocokan data dan kompleksitas model.
AIC dan BIC digunakan untuk mengevaluasi dan membandingkan model berdasarkan kompleksitas dan kecocokan; semakin kecil nilainya, semakin baik model tersebut.
Model regresi Poisson sering digunakan untuk memodelkan data yang berupa jumlah kejadian (count), seperti jumlah kecelakaan, jumlah kelahiran, atau jumlah pengunjung, dalam satuan waktu atau wilayah tertentu. Model ini sesuai diterapkan ketika variabel respon berupa bilangan bulat non-negatif dan diasumsikan mengikuti distribusi Poisson.
Distribusi Poisson memiliki bentuk fungsi probabilitas sebagai berikut:
\[ P(Y_i = y_i) = \frac{e^{-\lambda_i} \lambda_i^{y_i}}{y_i!} \]
di mana \(\lambda_i\) adalah rata-rata atau parameter intensitas dari distribusi Poisson untuk observasi ke-\(i\).
Dalam regresi Poisson, kita memodelkan hubungan antara nilai rata-rata \(\lambda_i\) dengan variabel prediktor \(x_i\) melalui fungsi link logaritmik:
\[ \log(\lambda_i) = \mathbf{x}_i^T \boldsymbol{\beta} \]
sehingga:
\[ \lambda_i = \exp(\mathbf{x}_i^T \boldsymbol{\beta}) \]
Model ini termasuk dalam keluarga Generalized Linear Model (GLM) dengan link fungsi log dan fungsi distribusi eksponensial dari Poisson.
Estimasi Parameter (MLE)
Untuk mengestimasi parameter \(\boldsymbol{\beta}\), digunakan pendekatan Maximum Likelihood Estimation (MLE). Fungsi log-likelihood untuk distribusi Poisson dapat ditulis sebagai:
\[ \ell(\boldsymbol{\beta}) = \sum_{i=1}^{n} \left[y_i \log(\lambda_i) - \lambda_i - \log(y_i!)\right] \]
dengan \(\lambda_i = \exp(\mathbf{x}_i^T \boldsymbol{\beta})\).
Estimasi dilakukan dengan metode iterasi (IRLS)
Estimasi MLE dilakukan secara iteratif menggunakan metode Iteratively Reweighted Least Squares (IRLS), yang merupakan pendekatan numerik untuk mencari solusi dari log-likelihood maksimum.
Langkah-langkah iterasi IRLS:
Tahap 1: Spesifikasi model
\[ \log(\lambda_i) = \mathbf{x}_i^T \boldsymbol{\beta} \Rightarrow \lambda_i = \exp(\mathbf{x}_i^T \boldsymbol{\beta}) \]
Tahap 2: Menyusun fungsi log-likelihood
\[ \ell(\boldsymbol{\beta}) = \sum_{i=1}^{n} \left[y_i \log(\lambda_i) - \lambda_i - \log(y_i!)\right] \]
Tahap 3: Formulasi iteratif
\[ \boldsymbol{\beta}^{(t+1)} = \left(\mathbf{X}^T \mathbf{W}^{(t)} \mathbf{X}\right)^{-1} \mathbf{X}^T \mathbf{W}^{(t)} \mathbf{z}^{(t)} \]
dengan:
Simulasi Data
set.seed(2025)
n <- 120
x1 <- rnorm(n, mean = 2, sd = 1)
x2 <- runif(n, min = -1, max = 1)
X <- cbind(1, x1, x2) # Menambahkan intercept
beta_true <- c(0.2, 0.6, -0.3)
eta <- X %*% beta_true
lambda <- exp(eta)
y <- rpois(n, lambda)
IRLS Manual Step-by-Step
# Inisialisasi
beta <- rep(0, ncol(X))
tol <- 1e-6
max_iter <- 100
for (i in 1:max_iter) {
eta <- X %*% beta
lambda <- exp(eta)
W <- diag(as.numeric(lambda))
z <- eta + (y - lambda) / lambda
beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z
if (sum(abs(beta_new - beta)) < tol) {
cat("Konvergen pada iterasi ke-", i, "\n")
break
}
beta <- beta_new
}
## Konvergen pada iterasi ke- 14
beta # hasil estimasi
## [,1]
## 0.05595847
## x1 0.65007233
## x2 -0.35209173
Perbandingan dengan glm()
model_glm <- glm(y ~ x1 + x2, family = poisson)
summary(model_glm)
##
## Call:
## glm(formula = y ~ x1 + x2, family = poisson)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.05596 0.11771 0.475 0.635
## x1 0.65007 0.04145 15.682 < 2e-16 ***
## x2 -0.35209 0.07883 -4.467 7.94e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 383.52 on 119 degrees of freedom
## Residual deviance: 125.17 on 117 degrees of freedom
## AIC: 493.88
##
## Number of Fisher Scoring iterations: 5
Metode IRLS berhasil mengestimasi parameter regresi Poisson dengan hasil yang sangat mendekati fungsi glm() di R. Estimasi parameter dari IRLS yaitu intercept = 0.0556, x1 = 0.6501, dan x2 = -0.3521, hampir identik dengan hasil glm(). IRLS konvergen dalam 14 iterasi, sedangkan glm() hanya butuh 5 iterasi karena lebih efisien. Kedua prediktor signifikan secara statistik (p < 0.001), dan nilai AIC sebesar 493.88. Ini menunjukkan bahwa IRLS bekerja dengan baik dan memberi pemahaman mendalam tentang proses estimasi di balik glm().
Pengujian hipotesis
Uji Wald
Untuk menguji \(H_0 = 0\)
coef_val <- coef(model)[3]
se_val <- summary(model)$coefficients[3, 2]
# Uji Wald
wald_z <- coef_val / se_val
wald_chisq <- wald_z^2
p_value <- 1 - pchisq(wald_chisq, df = 1)
cat("Z:", wald_z, "\nChi-Square:", wald_chisq, "\np-value:", p_value)
## Z: 4.820372
## Chi-Square: 23.23599
## p-value: 1.432904e-06
Uji Likelihood Ratio (Chi-Square)
model_null <- glm(y ~ x1, family = poisson, data = data)
anova(model_null, model, test = "Chisq")
## Analysis of Deviance Table
##
## Model 1: y ~ x1
## Model 2: y ~ x1 + x2
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 118 75.103
## 2 117 122.174 1 -47.071
Evaluasi Model
AIC(model)
## [1] 128.174
BIC(model)
## [1] 136.5365
Estimasi parameter regresi Poisson dilakukan dengan metode Maximum Likelihood Estimation (MLE). Untuk menguji signifikansi parameter, dilakukan Uji Wald dan Uji Likelihood Ratio. Nilai Z pada Uji Wald sebesar 4.820372 dengan nilai p-value sebesar 1.43e-06, menunjukkan bahwa parameter x2 signifikan secara statistik. Uji Likelihood Ratio antara model dengan prediktor x1 dan model lengkap (x1 dan x2) menghasilkan nilai deviance sebesar 75.103 dengan p-value signifikan, yang mengindikasikan bahwa penambahan x2 memperbaiki model secara signifikan. Evaluasi menggunakan AIC dan BIC menunjukkan nilai masing-masing 128.174 dan 136.5365. Nilai ini digunakan untuk membandingkan dan memilih model terbaik berdasarkan kompleksitas dan fit terhadap data.
Regresi logistik merupakan salah satu pendekatan dalam analisis statistik yang digunakan untuk menjelaskan hubungan antara suatu variabel dependen berskala biner (hanya memiliki dua kategori, misalnya: ya/tidak, sukses/gagal) dengan satu atau lebih variabel independen atau prediktor. Salah satu kelebihan dari regresi logistik adalah fleksibilitasnya dalam mengakomodasi berbagai jenis skala pengukuran pada variabel prediktor. Secara umum, skala pengukuran yang digunakan dalam regresi logistik dibedakan menjadi tiga, yaitu:
Skala Nominal: Variabel nominal adalah variabel yang kategorinya tidak memiliki urutan atau peringkat yang logis. Kategori dalam variabel ini hanya berfungsi sebagai pembeda tanpa menunjukkan tingkatan tertentu. Contoh yang umum adalah jenis kelamin (laki-laki dan perempuan), agama, atau warna favorit. Karena regresi logistik memerlukan data numerik, maka variabel nominal harus diubah menjadi bentuk numerik, biasanya melalui teknik dummy coding, yaitu membuat variabel biner (0 atau 1) untuk tiap kategori. Misalnya, untuk variabel jenis kelamin, kita dapat menetapkan perempuan = 1 dan laki-laki = 0.
Skala Ordinal: Variabel ordinal terdiri atas kategori-kategori yang memiliki urutan logis, tetapi belum tentu memiliki jarak yang sama antar kategori. Contohnya adalah tingkat pendidikan (SMA < Sarjana < Magister < Doktor), atau tingkat kepuasan (tidak puas < cukup puas < puas < sangat puas). Dalam praktiknya, variabel ordinal dapat diperlakukan dengan dua pendekatan berbeda:
Skala Rasio: Variabel rasio adalah variabel kuantitatif yang memiliki nol mutlak dan rasio antar nilai yang bermakna. Variabel ini memiliki semua sifat variabel interval, ditambah dengan adanya titik nol absolut yang menunjukkan tidak adanya atribut yang diukur. Contoh dari variabel rasio adalah pendapatan (misalnya dalam satuan juta rupiah per bulan), usia, atau jumlah anak. Karena sifatnya yang numerik dan memiliki skala yang lengkap, variabel ini dapat langsung digunakan dalam model regresi logistik tanpa perlu dilakukan transformasi khusus, kecuali jika distribusinya sangat tidak normal.
Seiring meningkatnya digitalisasi di sektor keuangan, penggunaan layanan online banking menjadi semakin umum. Namun, tidak semua individu memanfaatkan layanan ini secara merata. Oleh karena itu, dalam simulasi ini, kita ingin meneliti faktor-faktor yang mempengaruhi kemungkinan seseorang menggunakan layanan online banking, khususnya berdasarkan:
Saya membuat simulasi data sebanyak 750 observasi berdasarkan variabel-variabel tersebut.
# Jumlah observasi
n <- 750
# Simulasi variabel
job_status <- sample(c("Unemployed", "Employed", "Student", "Retired"),
n, replace = TRUE, prob = c(0.15, 0.5, 0.25, 0.1))
internet_use <- rpois(n, lambda = 7) # frekuensi penggunaan internet per minggu
age <- rnorm(n, mean = 35, sd = 10)
# Model logit
logit_p <- -1.5 + 0.6 * (job_status == "Employed") +
0.4 * (job_status == "Student") +
0.03 * internet_use - 0.02 * age
p <- 1 / (1 + exp(-logit_p))
set.seed(456)
online_banking <- rbinom(n, size = 1, prob = p)
# Gabungkan menjadi tibble
library(tibble)
## Warning: package 'tibble' was built under R version 4.3.3
sim_data <- tibble(online_banking, job_status, internet_use, age)
# Tampilkan sebagian data
head(sim_data)
## # A tibble: 6 × 4
## online_banking job_status internet_use age
## <int> <chr> <int> <dbl>
## 1 0 Unemployed 8 36.4
## 2 0 Employed 8 28.7
## 3 0 Student 9 39.0
## 4 1 Employed 6 26.3
## 5 0 Employed 15 43.7
## 6 0 Employed 6 14.2
Interpretasi Secara umum: Status bekerja dan frekuensi internet tinggi cenderung berasosiasi dengan penggunaan online banking.
Usia yang lebih muda juga tampak mendukung adopsi layanan ini.
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.3.3
##
## Attaching package: 'dplyr'
## The following object is masked from 'package:kableExtra':
##
## group_rows
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
sim_data %>%
group_by(online_banking) %>%
summarise(
Jumlah = n(),
Rata2_InternetUse = mean(internet_use),
Rata2_Age = mean(age)
)
## # A tibble: 2 × 4
## online_banking Jumlah Rata2_InternetUse Rata2_Age
## <int> <int> <dbl> <dbl>
## 1 0 594 6.93 34.7
## 2 1 156 7.18 34.2
Interpretasi Kelompok yang menggunakan online banking (online_banking = 1) berjumlah lebih banyak dibanding yang tidak menggunakan. Rata-rata usia mereka lebih muda dan memiliki frekuensi penggunaan internet lebih tinggi, yang mendukung hipotesis bahwa akses digital dan demografi muda mendorong adopsi layanan keuangan digital.
Sebuah lembaga perbankan ingin menganalisis faktor-faktor yang memengaruhi kemungkinan seseorang menggunakan layanan online banking. Salah satu variabel yang diteliti adalah status pekerjaan, yang terdiri dari empat kategori: Unemployed, Employed, Student, dan Retired.
Karena status pekerjaan bersifat kategorik (nominal), maka untuk kepentingan analisis regresi logistik, variabel ini diperlakukan sebagai variabel dummy. Kategori referensi yang digunakan adalah “Unemployed”.
Model logistik juga memasukkan variabel lain seperti frekuensi penggunaan internet per minggu (internet_use) dan usia (age) sebagai prediktor numerik.
Selanjutnya dilakukan regresi logistik untuk melihat pengaruh status pekerjaan (sebagai dummy) terhadap penggunaan online banking.
# Perlakuan status pekerjaan sebagai variabel dummy
sim_data_nominal <- sim_data %>%
mutate(
job_status = factor(job_status, levels = c("Unemployed", "Employed", "Student", "Retired"))
)
# Model regresi logistik
model_nominal <- glm(online_banking ~ job_status + internet_use + age,
data = sim_data_nominal, family = binomial)
summary(model_nominal)
##
## Call:
## glm(formula = online_banking ~ job_status + internet_use + age,
## family = binomial, data = sim_data_nominal)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.930209 0.477750 -4.040 5.34e-05 ***
## job_statusEmployed 0.809729 0.321117 2.522 0.0117 *
## job_statusStudent 0.501193 0.349101 1.436 0.1511
## job_statusRetired -0.266946 0.448853 -0.595 0.5520
## internet_use 0.039724 0.034901 1.138 0.2550
## age -0.006466 0.008791 -0.736 0.4620
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 766.94 on 749 degrees of freedom
## Residual deviance: 749.68 on 744 degrees of freedom
## AIC: 761.68
##
## Number of Fisher Scoring iterations: 4
Interpretasi Koefesien
Intercept = -1.9302 Rata-rata log odds untuk kelompok referensi (Unemployed) dengan nilai internet_use = 0 dan age = 0. Nilai ini tidak memiliki arti langsung secara praktis, karena usia dan penggunaan internet jarang 0, tapi tetap diperlukan dalam model.
job_statusEmployed = 0.8097, p = 0.0117 (*) Seseorang yang bekerja memiliki log odds menggunakan online banking lebih tinggi secara signifikan dibandingkan yang tidak bekerja (Unemployed). Odds Ratio = exp(0.8097) ≈ 2.25 → Lebih dari 2 kali lipat kemungkinan dibanding pengangguran.
job_statusStudent= 0.5012, p = 0.1511 (ns) Mahasiswa cenderung punya log odds lebih tinggi dibanding unemployed, tapi tidak signifikan secara statistik.
job_statusRetired Pensiunan memiliki kecenderungan lebih rendah dibanding unemployed, tidak signifikan.
internet_use = 0.0397, p = 0.2550 (ns) Frekuensi penggunaan internet tidak signifikan mempengaruhi penggunaan online banking dalam model ini.
age = -0.0065, p = 0.4620 (ns) Usia tidak berpengaruh signifikan terhadap peluang menggunakan layanan online banking.
Interpretasi Goodness of Fit
Null deviance: 766.94 Deviasi model tanpa prediktor
Residual deviance: 749.68 Deviasi model dengan semua prediktor
Penurunan deviance = 766.94 - 749.68 = 17.26 dengan df = 5
AIC = 761.68 Digunakan untuk perbandingan antar model (semakin kecil semakin baik).
Signifikansi Model
Hanya job_statusEmplyoed yang signifikan pada level 5%
Variabel lain tidak signikan, meskipun job_statusStudent hampir mendekati
Namun secara keseluruhan, model signifikan berdasarkan perbedaan deviance terhadap model null.
Kesimpulan Model menunjukkan bahwa:
Status pekerjaan sebagai “Employed” secara signifikan meningkatkan peluang seseorang menggunakan layanan online banking dibandingkan mereka yang tidak bekerja.
Variabel lain seperti status mahasiswa, pensiunan, usia, dan frekuensi internet belum menunjukkan pengaruh signifikan.
Model secara keseluruhan memberikan fit yang lebih baik daripada model tanpa prediktor, menunjukkan bahwa status pekerjaan memberi kontribusi informasi penting.
Dalam upaya memahami perilaku masyarakat terhadap layanan digital perbankan, salah satu faktor yang diduga berpengaruh adalah status pekerjaan. Status ini memiliki struktur yang dapat diasumsikan memiliki tingkatan, misalnya dari pengangguran (tidak bekerja) hingga bekerja penuh waktu atau pensiun. Oleh karena itu, pendekatan dengan memperlakukan variabel kategori ini sebagai skala numerik berjenjang (numeric rank) menjadi relevan.
Dengan asumsi bahwa terdapat urutan atau hirarki pada status pekerjaan, maka model regresi logistik dibentuk menggunakan perlakuan ordinal numerik terhadap variabel tersebut. Semakin tinggi nilai numeriknya, diinterpretasikan sebagai tingkatan pekerjaan yang lebih aktif atau mapan.
Pemodelan Regresi Logistik
library(dplyr)
# Ubah job_status menjadi numerik terurut
sim_data_numeric <- sim_data %>%
mutate(
job_status_numeric = case_when(
job_status == "Unemployed" ~ 1,
job_status == "Student" ~ 2,
job_status == "Retired" ~ 3,
job_status == "Employed" ~ 4
)
)
# Model regresi logistik dengan perlakuan ordinal
model_numeric <- glm(online_banking ~ job_status_numeric + internet_use + age,
data = sim_data_numeric, family = binomial)
summary(model_numeric)
##
## Call:
## glm(formula = online_banking ~ job_status_numeric + internet_use +
## age, family = binomial, data = sim_data_numeric)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.097605 0.464393 -4.517 6.28e-06 ***
## job_status_numeric 0.226756 0.084697 2.677 0.00742 **
## internet_use 0.034775 0.034759 1.000 0.31709
## age -0.005213 0.008782 -0.594 0.55279
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 766.94 on 749 degrees of freedom
## Residual deviance: 757.99 on 746 degrees of freedom
## AIC: 765.99
##
## Number of Fisher Scoring iterations: 4
Perbandingan Model
AIC_nominal = AIC(model_nominal); AIC_nominal
## [1] 761.6805
AIC_numeric = AIC(model_numeric); AIC_numeric
## [1] 765.9878
Interpretasi Koefesien
Intercept (-2.0976) mengindikasikan log odds dasar untuk individu pada tingkat status pekerjaan terendah (job_status_numeric = 1), dengan frekuensi penggunaan internet dan usia bernilai nol. Nilai ini tidak memiliki interpretasi praktis langsung, namun berfungsi sebagai titik awal perhitungan log odds.
job_status_numeric (0.2268, p = 0.0074) menunjukkan bahwa setiap kenaikan satu tingkat pada status pekerjaan (misalnya dari “Unemployed” ke “Student”) secara signifikan meningkatkan log odds seseorang menggunakan layanan online banking sebesar 0.2268. Dikonversi ke Odds Ratio, exp(0.2268) ≈ 1.254, yang berarti setiap kenaikan level status pekerjaan meningkatkan peluang penggunaan online banking sekitar 25.4%.
internet_use (0.0348, p = 0.3171) tidak menunjukkan pengaruh yang signifikan terhadap probabilitas penggunaan layanan online banking. Artinya, frekuensi penggunaan internet mingguan belum cukup kuat untuk memprediksi perilaku tersebut dalam model ini.
age (-0.0052, p = 0.5528) juga tidak berpengaruh signifikan. Penurunan log odds yang kecil ini tidak secara statistik berarti bahwa usia mempengaruhi penggunaan layanan online banking.
Interpretasi Goodness-of-Fit
Null deviance = 766.94 dan Residual deviance = 757.99 menunjukkan adanya penurunan sebesar 8.95, yang berarti model dengan prediktor memberikan informasi tambahan meskipun tidak terlalu besar.
AIC = 765.99, digunakan untuk menilai keseimbangan antara fit model dan kompleksitas. Dibandingkan dengan model sebelumnya (AIC = 761.68), model dengan perlakuan ordinal ini memiliki nilai AIC lebih tinggi, sehingga secara teknis model nominal lebih disukai berdasarkan AIC.
Signifikansi Model
Dalam model ini, hanya job_status_numeric yang memiliki pengaruh signifikan terhadap variabel dependen. Variabel lainnya tidak signifikan pada taraf 5%, sehingga peran prediktor selain status pekerjaan relatif lemah dalam menjelaskan variasi penggunaan layanan online banking.
Kesimpulan
Hasil analisis menunjukkan bahwa status pekerjaan yang dianggap memiliki tingkatan lebih tinggi berkorelasi positif dan signifikan terhadap peluang individu menggunakan layanan online banking. Meskipun penggunaan internet dan usia diikutsertakan dalam model, keduanya tidak terbukti signifikan. Bila dibandingkan dengan model sebelumnya, model dengan perlakuan nominal menunjukkan performa yang sedikit lebih baik berdasarkan nilai AIC, sehingga pendekatan kategorik dummy tampaknya lebih sesuai untuk kasus ini.
Goodness-of-Fit Model Untuk mengevaluasi kualitas model secara komprehensif, selain menggunakan deviance dan AIC, digunakan juga McFadden’s pseudo R², yang merupakan ukuran goodness-of-fit untuk model logistik. Nilai ini dihitung dengan membandingkan log-likelihood model penuh dengan model null (hanya intercept).
# Model null (tanpa prediktor)
nullmod <- glm(online_banking ~ 1, data = sim_data, family = binomial)
# McFadden's R²
r2_nominal <- 1 - (logLik(model_nominal) / logLik(nullmod))
r2_numeric <- 1 - (logLik(model_numeric) / logLik(nullmod))
list(
McFadden_R2_Nominal = r2_nominal,
McFadden_R2_Numeric = r2_numeric
)
## $McFadden_R2_Nominal
## 'log Lik.' 0.02250706 (df=6)
##
## $McFadden_R2_Numeric
## 'log Lik.' 0.01167538 (df=4)
Interpretasi
Nilai McFadden R² untuk model nominal sebesar 0.0225, sedangkan model ordinal sebesar 0.0117. Keduanya menunjukkan kemampuan jelas rendah dalam menjelaskan variabilitas data. Namun, model dengan perlakuan nominal memberikan kecocokan yang sedikit lebih baik dibanding model ordinal. Nilai ini mengindikasikan bahwa masih diperlukan variabel lain untuk meningkatkan performa model.
Visualisasi Prediksi (Model Nominal)
Prediksi probabilitas ini dihasilkan dari model regresi logistik dengan perlakuan status pekerjaan sebagai variabel dummy (nominal). Visualisasi ini menunjukkan hubungan antara frekuensi penggunaan internet (internet_use) dan prediksi probabilitas penggunaan layanan online banking, dikelompokkan berdasarkan status pekerjaan.
sim_data_nominal <- sim_data_nominal %>%
mutate(predicted = predict(model_nominal, type = "response"))
# Plot prediksi probabilitas berdasarkan status pekerjaan
ggplot(sim_data_nominal, aes(x = internet_use, y = predicted, color = job_status)) +
geom_point(alpha = 0.6) +
labs(
title = "Prediksi Probabilitas Penggunaan Online Banking (Model Nominal)",
x = "Frekuensi Penggunaan Internet per Minggu",
y = "Prediksi Probabilitas",
color = "Status Pekerjaan"
) +
theme_minimal() +
scale_color_manual(values = c(
"Unemployed" = "#E69F00", # orange
"Employed" = "#56B4E9", # blue
"Student" = "#009E73", # green
"Retired" = "#CC79A7" # pink
))
Visualisasi Prediksi Probabilitas (Model Ordinal/Numeric)
Visualisasi ini menampilkan prediksi probabilitas dari model regresi logistik dengan perlakuan ordinal (numeric rank) terhadap variabel status pekerjaan. Probabilitas dipetakan terhadap variabel pendapatan (income), dan diberi warna berdasarkan tingkat status pekerjaan yang telah dikodekan secara numerik.
library(ggplot2)
library(dplyr)
sim_data_numeric <- sim_data_numeric %>% mutate(predicted = predict(model_numeric, type = "response"))
# Plot probabilitas prediksi terhadap internet_use
ggplot(sim_data_numeric, aes(x = internet_use, y = predicted, color = as.factor(job_status_numeric))) +
geom_point(alpha = 0.6) +
labs(
title = "Prediksi Probabilitas Penggunaan Online Banking (Ordinal sebagai Numeric)",
x = "Frekuensi Penggunaan Internet per Minggu",
y = "Probabilitas Prediksi",
color = "Status Pekerjaan"
) +
theme_minimal() +
scale_color_manual(
values = c(
"1" = "#E41A1C", # merah terang
"2" = "#377EB8", # biru
"3" = "#4DAF4A", # hijau
"4" = "#984EA3" # ungu
),
labels = c("Unemployed", "Student", "Retired", "Employed")
)
Dalam analisis regresi logistik, proses pemilihan model merupakan langkah krusial untuk mendapatkan model yang paling sesuai dan efisien dalam memprediksi probabilitas kejadian suatu peristiwa biner. Model yang baik bukan hanya ditandai oleh ketepatan prediksi, tetapi juga oleh kesederhanaan, interpretabilitas, dan kesesuaian dengan teori atau tujuan penelitian.
Pemilihan model yang tepat juga dapat meningkatkan akurasi prediksi, menghindari overfitting, serta memperjelas kontribusi masing-masing prediktor terhadap respons. Untuk mencapai hal ini, biasanya digunakan pendekatan yang terstruktur dan metodologis, dengan mempertimbangkan baik aspek statistik maupun konteks substantif dari data yang dianalisis.
Secara umum, terdapat dua strategi dalam membangun model regresi logistik, yaitu pendekatan Confirmatory (konfirmatori) dan pendekatan Exploratory (eksploratori). Kedua pendekatan ini dapat digunakan secara terpisah maupun saling melengkapi, tergantung pada apakah peneliti telah memiliki teori atau masih dalam tahap eksplorasi hubungan antar variabel.
Contoh penggunaan: Misalnya, seorang peneliti ingin menguji apakah faktor X1 dan X2 berpengaruh terhadap keputusan pembelian produk. Dengan pendekatan konfirmatori, peneliti akan mengembangkan model berdasarkan teori perilaku konsumen, kemudian menguji apakah X1 dan X2 signifikan secara statistik.
Karateristik pendekatan eksploratori:
Proses pemilihan variabel dilakukan secara bertahap atau otomatis berdasarkan kriteria statistik tertentu.
Bertujuan untuk mengidentifikasi kombinasi variabel yang memberikan performa prediktif terbaik.
Pendekatan ini cocok digunakan untuk data besar atau studi awal eksplorasi fenomena.
Metode yang umum digunakan:
Forward Selection: Menambahkan variabel satu per satu ke dalam model berdasarkan kontribusi statistiknya.
Backward Elimination: Mengawali dengan semua variabel, kemudian mengeliminasi yang tidak signifikan.
Stepwise Selection: Kombinasi antara forward dan backward untuk memilih model terbaik.
Tujuan utama pendekatan ekspolatori
Menyusun model yang parsimonious (cukup sederhana) namun memiliki kinerja prediktif yang baik.
Cocok untuk studi eksplorasi, pengembangan model awal, atau ketika belum ada teori yang mapan.
Simulasi Data
# Memuat paket yang dibutuhkan
library(knitr)
library(dplyr)
library(ggplot2)
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
library(pROC)
## Warning: package 'pROC' was built under R version 4.3.3
## Type 'citation("pROC")' for a citation.
##
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
##
## cov, smooth, var
library(DescTools)
## Warning: package 'DescTools' was built under R version 4.3.3
set.seed(456)
# Jumlah sampel
n <- 250
# Prediktor X1 dan X3 bertipe numerik kontinu
x1 <- rnorm(n, mean = 2, sd = 1)
x3 <- rnorm(n, mean = 0, sd = 2)
# Prediktor X2 bertipe biner (simulasi variabel dummy)
x2 <- sample(0:1, n, replace = TRUE, prob = c(0.4, 0.6))
# Menyusun persamaan linear prediktor
lin_pred <- -1 + 0.7 * x1 + 0.4 * x2 - 0.6 * x3
# Menghitung probabilitas menggunakan fungsi logistik
p <- plogis(lin_pred)
# Menghasilkan variabel respon biner berdasarkan probabilitas
y <- rbinom(n, size = 1, prob = p)
# Menyusun data frame
df <- data.frame(
y = factor(y),
x1 = x1,
x2 = factor(x2),
x3 = x3
)
# Tampilkan sebagian data
head(df)
## y x1 x2 x3
## 1 1 0.6564786 1 -0.4998270
## 2 1 2.6217756 1 0.1426120
## 3 1 2.8008747 1 0.5001287
## 4 1 0.6111076 0 -2.4578501
## 5 1 1.2856431 1 -0.1172847
## 6 1 1.6759389 1 -0.7275040
Pemodelan regresi logistik penuh
# Pemodelan regresi logistik penuh
model_lengkap <- glm(y ~ x1 + x2 + x3, data = df, family = binomial)
# Menampilkan ringkasan model
summary(model_lengkap)
##
## Call:
## glm(formula = y ~ x1 + x2 + x3, family = binomial, data = df)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.81302 0.37873 -2.147 0.031818 *
## x1 0.57650 0.15585 3.699 0.000216 ***
## x21 0.46931 0.29119 1.612 0.107023
## x3 -0.48096 0.08612 -5.585 2.34e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 337.30 on 249 degrees of freedom
## Residual deviance: 285.31 on 246 degrees of freedom
## AIC: 293.31
##
## Number of Fisher Scoring iterations: 4
Untuk menyusun model regresi logistik yang optimal, teknik stepwise selection dapat digunakan untuk memilih kombinasi variabel prediktor terbaik secara otomatis. Pendekatan ini dilakukan dengan menambahkan (forward), menghapus (backward), atau mengombinasikan keduanya (both) berdasarkan kriteria statistik seperti AIC (Akaike Information Criterion).
# Model nol (hanya intercept, tanpa prediktor)
model_nol <- glm(y ~ 1, data = df, family = binomial)
# Seleksi maju (forward selection)
maju <- step(model_nol, direction = "forward", scope = formula(model_lengkap), trace = FALSE)
# Seleksi mundur (backward elimination)
mundur <- step(model_lengkap, direction = "backward", trace = FALSE)
# Seleksi dua arah (forward + backward / stepwise)
gabungan <- step(model_nol, direction = "both", scope = formula(model_lengkap), trace = FALSE)
Perbandingan AIC antar Metode Untuk membandingkan kebaikan model hasil seleksi, digunakan nilai AIC. Model dengan AIC lebih kecil dianggap lebih baik karena memberikan keseimbangan antara kompleksitas dan kualitas prediksi.
AIC(model_lengkap, maju, mundur, gabungan)
## df AIC
## model_lengkap 4 293.3136
## maju 4 293.3136
## mundur 4 293.3136
## gabungan 4 293.3136
Untuk menilai kemampuan model dalam membedakan antara kelas 1 dan 0, digunakan kurva ROC (Receiver Operating Characteristic). ROC menggambarkan hubungan antara Sensitivity (True Positive Rate) dan 1 − Specificity (False Positive Rate). Area di bawah kurva (AUC) menjadi ukuran numerik performa model klasifikasi.
# Prediksi probabilitas dari model terbaik (hasil stepwise)
prediksi_probabilitas <- predict(gabungan, type = "response")
# Membuat objek ROC
roc_kurva <- roc(df$y, prediksi_probabilitas)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
# Visualisasi kurva ROC
plot(roc_kurva,
main = "Kurva ROC",
col = "darkorange",
lwd = 2)
# Tambahkan garis diagonal sebagai referensi
abline(a = 0, b = 1, lty = 2, col = "gray")
# Tambahkan legend
legend("bottomright",
legend = c("Model Stepwise"),
col = c("darkorange"),
lwd = 2)
auc(roc_kurva)
## Area under the curve: 0.7539
Interpretasi
Nilai AUC (Area Under the Curve) sebesar 0.7539 menunjukkan bahwa model memiliki kemampuan yang cukup baik dalam membedakan antara dua kelas (respon 0 dan 1). Secara umum, nilai AUC di atas 0.7 dianggap memadai, yang berarti bahwa dalam sekitar 75% kasus, model dapat membedakan dengan benar antara observasi yang termasuk dalam kelas positif dan kelas negatif. Dengan demikian, model ini dapat dikatakan memiliki kinerja klasifikasi yang layak dan dapat digunakan sebagai dasar untuk pengambilan keputusan prediktif, meskipun masih terdapat ruang untuk perbaikan lebih lanjut.
Sebagai alternatif dari R-squared pada regresi linear, regresi logistik menggunakan pseudo R-squared untuk mengukur seberapa baik model menjelaskan variasi dalam data. Beberapa ukuran yang umum digunakan adalah Cox & Snell, Nagelkerke, dan McFadden
PseudoR2(gabungan, which = c("CoxSnell", "Nagelkerke", "McFadden"))
## CoxSnell Nagelkerke McFadden
## 0.1877492 0.2535258 0.1541254
Nilai pseudo R-squared ini mengindikasikan bahwa model dapat menjelaskan sekitar 15%–25% variasi dalam respons, yang dianggap cukup baik untuk model klasifikasi dengan data sosial atau observasional.
Untuk menilai performa prediksi, dibuat tabel klasifikasi (confusion matrix) dengan cutoff probabilitas 0.5. Evaluasi mencakup akurasi, sensitivitas, spesifisitas, serta metrik tambahan lainnya.
# Mengklasifikasikan berdasarkan probabilitas prediksi
kelas_prediksi <- ifelse(prediksi_probabilitas >= 0.5, 1, 0)
# Membuat tabel perbandingan antara prediksi dan data aktual
conf_matrix <- table(Prediksi = kelas_prediksi, Aktual = df$y)
print(conf_matrix)
## Aktual
## Prediksi 0 1
## 0 53 23
## 1 48 126
# Ambil nilai-nilai dari confusion matrix
TN <- conf_matrix[1,1]
FP <- conf_matrix[2,1]
FN <- conf_matrix[1,2]
TP <- conf_matrix[2,2]
# Hitung metrik
accuracy <- (TP + TN) / sum(conf_matrix)
sensitivity <- TP / (TP + FN) # Recall
specificity <- TN / (TN + FP)
precision <- TP / (TP + FP)
f1_score <- 2 * (precision * sensitivity) / (precision + sensitivity)
# Tampilkan hasil
cat("Akurasi :", round(accuracy, 3), "\n")
## Akurasi : 0.716
cat("Sensitivitas:", round(sensitivity, 3), "\n")
## Sensitivitas: 0.846
cat("Spesifisitas:", round(specificity, 3), "\n")
## Spesifisitas: 0.525
cat("Presisi :", round(precision, 3), "\n")
## Presisi : 0.724
cat("F1-Score :", round(f1_score, 3), "\n")
## F1-Score : 0.78
Model klasifikasi menunjukkan kinerja yang cukup baik dengan akurasi sebesar 71,6%, menandakan bahwa sebagian besar prediksi sudah tepat. Nilai sensitivitas yang tinggi (84,6%) menunjukkan kemampuan model yang kuat dalam mendeteksi kasus positif, sementara presisi sebesar 72,4% mengindikasikan bahwa sebagian besar prediksi positif memang benar. Meskipun spesifisitasnya masih rendah (52,5%), secara keseluruhan model cukup andal, terutama jika fokus utama adalah mendeteksi kelas positif.
Tujuan
Bagian ini bertujuan untuk menunjukkan bagaimana membandingkan beberapa model regresi logistik dengan menggunakan ukuran statistik yang umum, yaitu: - AIC (Akaike Information Criterion) - Deviance (−2 log-likelihood) - Likelihood-Ratio Test Selain itu, dijelaskan pula pentingnya prinsip parsimony, yakni memilih model yang cukup sederhana namun tetap mampu menjelaskan data secara efektif.
Simulasi Data
library(MASS)
library(broom)
## Warning: package 'broom' was built under R version 4.3.3
library(DescTools)
set.seed(456)
n <- 500
z1 <- rnorm(n)
z2 <- rbinom(n, 1, 0.5)
z3 <- rnorm(n)
linear_pred <- -1.3 + 1.8 * z1 - 0.1 * z2 + 0.75 * z3
p <- 1 / (1 + exp(-linear_pred))
target <- rbinom(n, 1, p)
data_log <- data.frame(y = as.factor(target), z1, z2, z3)
Pembuatan Model Dibuat tiga model regresi logistik dengan kompleksitas bertingkat, dari model sederhana hingga model lengkap:
model_a <- glm(y ~ z1, data = data_log, family = binomial)
model_b <- glm(y ~ z1 + z2, data = data_log, family = binomial)
model_c <- glm(y ~ z1 + z2 + z3, data = data_log, family = binomial)
Perbandingan AIC dan Deviance Untuk memilih model terbaik, dilakukan perbandingan berdasarkan nilai AIC dan deviance:
model_band <- data.frame(
Model = c("Model A", "Model B", "Model C"),
AIC = c(AIC(model_a), AIC(model_b), AIC(model_c)),
Deviance = c(deviance(model_a), deviance(model_b), deviance(model_c))
)
model_band
## Model AIC Deviance
## 1 Model A 463.0457 459.0457
## 2 Model B 464.7974 458.7974
## 3 Model C 431.5352 423.5352
anova(model_a, model_b, test = "LRT")
## Analysis of Deviance Table
##
## Model 1: y ~ z1
## Model 2: y ~ z1 + z2
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 498 459.05
## 2 497 458.80 1 0.24835 0.6182
anova(model_b, model_c, test = "LRT")
## Analysis of Deviance Table
##
## Model 1: y ~ z1 + z2
## Model 2: y ~ z1 + z2 + z3
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 497 458.80
## 2 496 423.54 1 35.262 2.882e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Dalam pemodelan regresi logistik, kompleksitas model sering kali dikaitkan dengan peningkatan kecocokan terhadap data. Namun, model yang terlalu kompleks berisiko menimbulkan overfitting dan dapat menyulitkan interpretasi. Oleh karena itu, prinsip parsimony menekankan pentingnya memilih model yang relatif sederhana namun tetap mampu menggambarkan struktur data secara memadai.
Secara umum, jika penambahan parameter tidak memberikan penurunan nilai AIC yang signifikan, maka model yang lebih sederhana sebaiknya dipilih. Prinsip ini mendukung penggunaan model yang efisien dan tidak berlebihan dalam jumlah parameter.
Rumus dan Penjelasan
Rumus Deviance \[ D = -2 \left[ \log L(\text{model}) - \log L(\text{model saturasi}) \right] \] Deviance mengukur sejauh mana suatu model berbeda dari model saturasi (model dengan kecocokan sempurna terhadap data). Semakin kecil nilai deviance, maka model tersebut dianggap semakin mendekati data aktual. Dengan demikian, deviance yang rendah mencerminkan prediksi model yang lebih akurat.
Rumus Likelihood-Ratio Test \[ G^2 = -2(\log L_0 - \log L_1) \] Statistik Likelihood Ratio (\(G^2\)) digunakan untuk mengevaluasi apakah penambahan variabel dalam suatu model memberikan peningkatan kecocokan model secara signifikan. Nilai \(G^2\) yang besar dan diikuti oleh nilai p-value yang kecil menunjukkan bahwa model yang lebih kompleks secara statistik lebih unggul dibandingkan model yang lebih sederhana.
Setelah model regresi logistik dibangun, langkah penting berikutnya adalah mengevaluasi performa prediktifnya. Salah satu pendekatan yang umum digunakan adalah membentuk tabel klasifikasi (confusion matrix) berdasarkan prediksi model dan membandingkannya dengan data aktual.
Klasifikasi dilakukan dengan menetapkan ambang batas (cut-off) probabilitas, misalnya 0.5. Observasi dengan probabilitas di atas cut-off diklasifikasikan sebagai kelas positif, dan sisanya sebagai kelas negatif.
Berikut ini adalah implementasi evaluasi tabel klasifikasi
tanpa menggunakan paket tambahan seperti
caret
, melainkan dengan perhitungan manual.
# Menghasilkan probabilitas prediksi dari model logistik
pred_prob <- predict(model_c, type = "response")
# Klasifikasi berdasarkan cutoff 0.5
pred_class <- ifelse(pred_prob >= 0.5, 1, 0)
# Data aktual
actual_class <- as.numeric(as.character(data_log$y)) # pastikan dalam bentuk numerik
# Matriks klasifikasi
conf_matrix <- matrix(0, nrow = 2, ncol = 2,
dimnames = list("Predicted" = c("0", "1"),
"Actual" = c("0", "1")))
# Mengisi confusion matrix
for (i in 1:length(actual_class)) {
predicted <- pred_class[i]
actual <- actual_class[i]
conf_matrix[predicted + 1, actual + 1] <- conf_matrix[predicted + 1, actual + 1] + 1
}
# Menampilkan confusion matrix
conf_matrix
## Actual
## Predicted 0 1
## 0 293 61
## 1 36 110
# Ekstrak nilai dari confusion matrix
TN <- conf_matrix["0", "0"]
FP <- conf_matrix["1", "0"]
FN <- conf_matrix["0", "1"]
TP <- conf_matrix["1", "1"]
# Hitung metrik performa
akurasi <- (TP + TN) / sum(conf_matrix)
sensitivitas <- TP / (TP + FN) # Recall
spesifisitas <- TN / (TN + FP)
presisi <- TP / (TP + FP)
f1_score <- 2 * (presisi * sensitivitas) / (presisi + sensitivitas)
# Tampilkan hasil evaluasi
cat("Akurasi :", round(akurasi, 3), "\n")
## Akurasi : 0.806
cat("Sensitivitas:", round(sensitivitas, 3), "\n")
## Sensitivitas: 0.643
cat("Spesifisitas:", round(spesifisitas, 3), "\n")
## Spesifisitas: 0.891
cat("Presisi :", round(presisi, 3), "\n")
## Presisi : 0.753
cat("F1-Score :", round(f1_score, 3), "\n")
## F1-Score : 0.694
Interpretasi
Berdasarkan hasil evaluasi model klasifikasi, diperoleh nilai akurasi sebesar 0.806, yang menunjukkan bahwa sekitar 80.6% dari seluruh prediksi model telah sesuai dengan data aktual. Sensitivitas sebesar 0.643 mengindikasikan bahwa model mampu mengidentifikasi 64.3% dari kasus positif dengan benar, sedangkan spesifisitas sebesar 0.891 menandakan kemampuan yang sangat baik dalam mengenali kasus negatif. Presisi (0.753) mengukur bahwa dari seluruh prediksi positif, 75.3% adalah benar-benar positif. Adapun nilai F1-Score sebesar 0.694 merepresentasikan harmoni antara presisi dan sensitivitas, yang penting dalam konteks klasifikasi tidak seimbang. Secara keseluruhan, model memiliki kinerja yang cukup baik, terutama dalam mendeteksi kelas negatif, meskipun performa terhadap kelas positif masih dapat ditingkatkan.
Sensitivitas dan Spesifisitas Untuk mengevaluasi secara lebih spesifik kemampuan model dalam klasifikasi biner, dua metrik utama yang perlu diperhatikan adalah sensitivitas dan spesifisitas. Kedua ukuran ini merepresentasikan seberapa baik model mengenali kelas positif dan negatif secara terpisah.
Sensitivitas Menggambarkan kelas positif yang berhasil teridentifikasi dengan benar oleh model. Ukuran ini dikenal juga sebagai True Positive Rate dan dirumuskan sebagai:
\[ \text{Sensitivity} = \frac{\text{TP}}{\text{TP} + \text{FN}} \]
Spesifisitas Menggambarkan kelas negatif yang dikenali secara tepat. Ini disebut pula sebagai True Negative Rate, dengan rumus sebagai berikut:
\[ \text{Specificity} = \frac{\text{TN}}{\text{TN} + \text{FP}} \]
Kedua metrik tersebut sangat penting untuk memahami ketidakseimbangan kinerja model dalam mendeteksi dua kategori kelas, terutama ketika fokus utama berada pada deteksi kejadian positif yang mungkin jarang terjadi namun kritis untuk diidentifikasi secara akurat.
Perhitungan sensitivitas dan spesifisitas dapat dilakukan secara manual maupun menggunakan kode berikut:
# Perhitungan manual telah dilakukan sebelumnya
cat("Sensitivitas:", round(sensitivitas, 3), "\n")
## Sensitivitas: 0.643
cat("Spesifisitas:", round(spesifisitas, 3), "\n")
## Spesifisitas: 0.891
Kesimpulan
Deviance yang rendah mencerminkan bahwa model memiliki tingkat kesalahan prediksi yang minim dibandingkan model referensi yang sempurna.
AIC digunakan untuk menilai keseimbangan antara ketepatan model dan jumlah parameter yang digunakan, sehingga membantu menghindari model yang terlalu kompleks.
Uji Likelihood Ratio berperan dalam menentukan apakah penambahan variabel pada model memberikan peningkatan signifikan secara statistik.
Tabel klasifikasi memberikan gambaran kuantitatif mengenai akurasi model dalam membedakan antara kategori aktual dan hasil prediksi.
Prinsip Parsimony mengarahkan pemilihan model ke arah struktur yang lebih sederhana apabila tidak terdapat perbedaan kinerja yang substansial dengan model yang lebih rumit.
Kurva ROC (Receiver Operating Characteristic) merupakan representasi grafis yang digunakan dalam evaluasi kinerja model klasifikasi biner. Kurva ini menunjukkan hubungan antara True Positive Rate (Sensitivitas) dan False Positive Rate (1 − Spesifisitas) pada berbagai nilai ambang klasifikasi. Semakin tinggi kemampuan model dalam membedakan antara kelas positif dan negatif, maka kurva ROC akan semakin menjauhi garis diagonal dan mendekati sudut kiri atas grafik.
Sumbu Y (Sensitivitas): Kelas positif yang berhasil dikenali model secara benar:
\[ \text{Sensitivitas} = \frac{TP}{TP + FN} \]
Sumbu X (1 − Spesifisitas): Proporsi kesalahan klasifikasi terhadap kelas negatif (False Positive Rate):
\[ \text{FPR} = \frac{FP}{FP + TN} = 1 - \text{Spesifisitas} \]
Garis diagonal (y = x) menunjukkan performa acak (random classifier), yang berarti model tidak memiliki kemampuan diskriminatif.
Kurva yang mendekati sudut kiri atas merefleksikan performa klasifikasi yang baik karena menunjukkan kombinasi sensitivitas tinggi dan false positive rate yang rendah.
Perubahan nilai ambang klasifikasi (cut-off) memengaruhi bentuk kurva ROC:
Model yang ideal secara teoritis akan menghasilkan kurva ROC dengan ciri:
AUC juga dikenal sebagai ukuran concordance, yaitu probabilitas bahwa model memberikan skor yang lebih tinggi untuk kasus positif dibanding kasus negatif secara acak.
library(pROC)
set.seed(321)
# Simulasi variabel prediktor
x1 <- rnorm(300, mean = 0.5, sd = 1.3)
x2 <- rbinom(300, size = 1, prob = 0.6)
x3 <- rnorm(300, mean = -1, sd = 1)
# Model linier dengan parameter berbeda
lin_mod <- 0.8 - 1.4 * x1 + 1.2 * x2 + 0.5 * x3
p <- 1 / (1 + exp(-lin_mod))
y <- rbinom(300, 1, prob = p)
# Data frame
data <- data.frame(y = factor(y), x1, x2, x3)
# Model logistik
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
# Prediksi probabilitas
pred <- predict(model, type = "response")
# Kurva ROC
roc_curve <- roc(data$y, pred)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
# Plot ROC dengan warna unik dan legend
plot(roc_curve,
main = "Kurva ROC - Model Simulasi",
col = "#4B0082", # warna ungu gelap (indigo)
lwd = 2)
abline(a = 0, b = 1, col = "gray", lty = 2)
legend("bottomright",
legend = c("Model Logistik Simulasi"),
col = c("#4B0082"),
lwd = 2)
auc(roc_curve)
## Area under the curve: 0.8614
# Daftar nilai threshold yang akan diuji
thresholds <- seq(0.05, 0.95, by = 0.05)
# Inisialisasi data frame untuk menyimpan hasil
results <- data.frame(Threshold = thresholds)
# Hitung Sensitivitas untuk setiap threshold
results$Sensitivity <- sapply(thresholds, function(t) {
pred_class <- ifelse(pred >= t, 1, 0)
cm <- table(Pred = pred_class, Obs = as.numeric(as.character(data$y)))
TP <- cm["1", "1"]
FN <- cm["0", "1"]
TP / (TP + FN)
})
# Hitung Spesifisitas untuk setiap threshold
results$Specificity <- sapply(thresholds, function(t) {
pred_class <- ifelse(pred >= t, 1, 0)
cm <- table(Pred = pred_class, Obs = as.numeric(as.character(data$y)))
TN <- cm["0", "0"]
FP <- cm["1", "0"]
TN / (TN + FP)
})
# Cetak hasil evaluasi sensitivitas dan spesifisitas
print(results)
## Threshold Sensitivity Specificity
## 1 0.05 1.0000000 0.1313869
## 2 0.10 0.9938650 0.2116788
## 3 0.15 0.9631902 0.3138686
## 4 0.20 0.9570552 0.3795620
## 5 0.25 0.9386503 0.4598540
## 6 0.30 0.9325153 0.5255474
## 7 0.35 0.9079755 0.5912409
## 8 0.40 0.8773006 0.6569343
## 9 0.45 0.8588957 0.7153285
## 10 0.50 0.8282209 0.7737226
## 11 0.55 0.7730061 0.7956204
## 12 0.60 0.7177914 0.8394161
## 13 0.65 0.6871166 0.8759124
## 14 0.70 0.6257669 0.9051095
## 15 0.75 0.5705521 0.9343066
## 16 0.80 0.5092025 0.9416058
## 17 0.85 0.4294479 0.9416058
## 18 0.90 0.2699387 0.9635036
## 19 0.95 0.1104294 0.9854015
Penjelasan Umum Kurva Precision-Recall (PR) merupakan salah satu metode evaluasi performa model klasifikasi, yang sangat direkomendasikan terutama saat kita menghadapi data dengan ketidakseimbangan kelas (class imbalance), di mana jumlah observasi pada kelas positif jauh lebih sedikit dibanding kelas negatif.
\[ \text{Precision} = \frac{TP}{TP + FP} \]
\[ \text{Recall} = \frac{TP}{TP + FN} \] Keterangan Notasi:
\(TP\): True Positive (prediksi positif yang benar)
\(FP\): False Positive (prediksi positif yang salah)
\(FN\): False Negative (prediksi negatif yang salah)
Catatan: PR Curve lebih representatif dibanding ROC Curve untuk kasus data yang tidak seimbang karena tidak memperhitungkan true negatives.
\[ \text{Baseline AUPRC} = \frac{\text{Jumlah Kelas Positif}}{\text{Total Observasi}} \] AUPRC sangat bermanfaat untuk menilai model dalam kasus ketidakseimbangan kelas, karena fokus pada kelas positif.
Aspek | ROC Curve | Precision-Recall Curve |
---|---|---|
Fokus | Semua kelas | Kelas positif saja |
Kuat di | Data seimbang | Data tidak seimbang |
Sumbu Y | Sensitivitas (Recall) | Presisi |
Sumbu X | 1 - Spesifisitas (FPR) | Recall |
# Muat paket yang diperlukan
library(PRROC)
## Warning: package 'PRROC' was built under R version 4.3.3
## Loading required package: rlang
## Warning: package 'rlang' was built under R version 4.3.3
library(rlang)
# Inisialisasi seed untuk replikasi
set.seed(2025)
# Simulasi tiga prediktor dengan distribusi normal
x1 <- rnorm(300)
x2 <- rnorm(300, 2, 0.4)
x3 <- rnorm(300)
# Rumus logit dan prediksi probabilitas
lin_pred <- -1 + 0.9 * x1 + 0.6 * x2 + 0.5 * x3
p <- 1 / (1 + exp(-lin_pred))
# Hasil klasifikasi biner
y <- rbinom(300, 1, p)
# Gabungkan data menjadi satu kerangka data
data <- data.frame(y = y, x1 = x1, x2 = x2, x3 = x3)
# Bangun model regresi logistik
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
# Prediksi probabilitas
prob <- predict(model, type = "response")
# Buat objek PR Curve
pr <- pr.curve(scores.class0 = prob[data$y == 1],
scores.class1 = prob[data$y == 0],
curve = TRUE)
# Visualisasikan
plot(pr)
Dalam regresi logistik, koefisien determinasi (\(R^2\)) seperti yang dikenal pada regresi linear tidak dapat digunakan secara langsung karena bentuk fungsi likelihood yang berbeda. Oleh karena itu, dikembangkanlah ukuran alternatif yang disebut pseudo R-squared.
Tujuan dari penggunaan pseudi R-squareq
Memberikan indikasi relatif seberapa baik model logistik menjelaskan variabilitas data biner yang diamati.
Membandingkan kebaikan model antara model penuh dengan model null (tanpa prediktor).
Menilai apakah penambahan prediktor meningkatkan kecocokan model terhadap data.
Memberi bantuan kuantitatif dalam pemilihan model (selain menggunakan AIC/BIC atau deviance).
- $R^2_{\text{CoxSnell}}$
- $R^2_{\text{McFadden}}$
Simulasi Data
# Set seed agar hasil dapat direproduksi
set.seed(58)
# Jumlah observasi
n <- 500
# Simulasi prediktor x1, x2, dan x3
x1 <- rnorm(n)
x2 <- rbinom(n, 1, 0.4)
x3 <- rnorm(n)
# Persamaan linier logit dan probabilitas
lin_pred <- -0.8 + 1.5 * x1 - 0.3 * x2 + 0.9 * x3
p <- 1 / (1 + exp(-lin_pred))
# Respon biner y
y <- rbinom(n, 1, p)
# Gabungkan ke dalam dataframe
data <- data.frame(y = as.factor(y), x1 = x1, x2 = x2, x3 = x3)
Model Logistik dan Null Model
# Model penuh (menggunakan semua prediktor)
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
model
##
## Call: glm(formula = y ~ x1 + x2 + x3, family = binomial, data = data)
##
## Coefficients:
## (Intercept) x1 x2 x3
## -0.8641 1.6928 -0.4486 0.8968
##
## Degrees of Freedom: 499 Total (i.e. Null); 496 Residual
## Null Deviance: 639.7
## Residual Deviance: 442.4 AIC: 450.4
# Model null (hanya intercept, tanpa prediktor)
model_null <- glm(y ~ 1, data = data, family = binomial)
model_null
##
## Call: glm(formula = y ~ 1, family = binomial, data = data)
##
## Coefficients:
## (Intercept)
## -0.6722
##
## Degrees of Freedom: 499 Total (i.e. Null); 499 Residual
## Null Deviance: 639.7
## Residual Deviance: 639.7 AIC: 641.7
Interpretasi
Berdasarkan hasil regresi logistik, model penuh yang menggunakan
prediktor x1
, x2
, dan x3
menunjukkan bahwa x1
dan x3
memiliki
kontribusi positif terhadap probabilitas kejadian (dengan koefisien
masing-masing 1.6928 dan 0.8968), sedangkan x2
memiliki
kontribusi negatif (-0.4486). Perbedaan antara null deviance (639.7) dan
residual deviance (442.4) menunjukkan bahwa model dengan prediktor
memberikan peningkatan kecocokan yang signifikan dibanding model tanpa
prediktor. Nilai \(R^2_{\text{CoxSnell}}\) sebesar 0.2698 dan
\(R^2_{\text{McFadden}}\) sebesar
0.2586 menunjukkan bahwa model menjelaskan sebagian variasi dalam data
dengan kualitas yang tergolong cukup baik untuk model
logistik biner.
Likelihood dan Rumus Ukuran Pseudo \(R^2\) dihitung berdasarkan nilai likelihood sebagai berikut:
Rumus Cox & Snell:
\[ R^2_{\text{CoxSnell}} = 1 - \left( \frac{L_0}{L_M} \right)^{\frac{2}{n}} \]
Rumus McFadden:
\[ R^2_{\text{McFadden}} = 1 - \frac{\log L_M}{\log L_0} \] Dengan:
Perhitungan Manual dalam R
# Ambil nilai log-likelihood
logL0 <- logLik(model_null)
logLM <- logLik(model)
# Konversi ke skala likelihood
L0 <- exp(logL0)
LM <- exp(logLM)
# Ambil jumlah data
n <- nobs(model)
# Hitung pseudo R-squared
cox_snell <- 1 - (L0 / LM)^(2 / n)
mcfadden <- 1 - (as.numeric(logLM) / as.numeric(logL0))
# Tampilkan dalam data.frame
r2 <- data.frame(
R2_CoxSnell = cox_snell,
R2_McFadden = mcfadden
)
r2
## R2_CoxSnell R2_McFadden
## 1 0.3260752 0.3084545
Perhitungan Otomatis dengan beberapa Package
library(pscl)
## Warning: package 'pscl' was built under R version 4.3.3
## Classes and Methods for R originally developed in the
## Political Science Computational Laboratory
## Department of Political Science
## Stanford University (2002-2015),
## by and under the direction of Simon Jackman.
## hurdle and zeroinfl functions by Achim Zeileis.
pR2(model)
## fitting null model for pseudo-r2
## llh llhNull G2 McFadden r2ML r2CU
## -221.1908043 -319.8499841 197.3183597 0.3084545 0.3260752 0.4517554
library(DescTools)
PseudoR2(model, which = "all")
## McFadden McFaddenAdj CoxSnell Nagelkerke AldrichNelson
## 0.3084545 0.2959487 0.3260752 0.4517554 0.2829674
## VeallZimmermann Efron McKelveyZavoina Tjur AIC
## 0.5041394 0.3424673 0.5002325 0.3508718 450.3816086
## BIC logLik logLik0 G2
## 467.2400410 -221.1908043 -319.8499841 197.3183597
Interpretasi Hasil evaluasi model logistik dengan
menggunakan paket pscl
dan DescTools
menunjukkan bahwa nilai pseudo R-squared model cukup bervariasi
tergantung metode yang digunakan. Nilai \(R^2_{\text{McFadden}}\) sebesar 0.3085
menunjukkan bahwa model menjelaskan sekitar 30.85% dari informasi
log-likelihood dibandingkan model null. Sementara itu, nilai \(R^2_{\text{CoxSnell}}\) sebesar 0.3261 dan
\(R^2_{\text{Nagelkerke}}\) sebesar
0.4518 memberikan gambaran bahwa model memiliki kualitas prediksi yang
sedang hingga baik. Nilai \(R^2\)
lainnya seperti Tjur (0.5041) dan McKelvey-Zavoina (0.5002) menunjukkan
bahwa setidaknya setengah dari variasi dalam probabilitas kejadian dapat
dijelaskan oleh model. Secara keseluruhan, metrik-metrik ini
mengindikasikan bahwa model memiliki kecocokan yang cukup memadai
terhadap data.
Distribusi multinomial merupakan perluasan dari distribusi binomial ke dalam konteks yang memiliki lebih dari dua kategori hasil. Misalkan terdapat \(X_1, X_2, \dots, X_k\) sebagai banyaknya observasi untuk masing-masing dari \(k\) kategori, maka probabilitas gabungan dari distribusi multinomial adalah:
\[ P(X_1 = x_1, \dots, X_k = x_k) = \frac{n!}{x_1! x_2! \dots x_k!} \cdot p_1^{x_1} p_2^{x_2} \dots p_k^{x_k} \]
dengan: - \(\sum_{i=1}^{k} x_i = n\) - \(\sum_{i=1}^{k} p_i = 1\)
Sebuah penelitian dilakukan terhadap 15 responden yang diminta untuk memilih satu dari empat jenis minuman favorit: Teh (T), Kopi (K), Jus (J), dan Soda (S). Hasil observasi yang diperoleh adalah sebagai berikut:
Asumsikan bahwa probabilitas teoritis preferensi minuman di populasi adalah sebagai berikut:
Pertanyaan:
Berapa peluang bahwa dalam 15 responden, sebanyak 6 orang memilih Teh, 4
orang memilih Kopi, 3 orang memilih Jus, dan 2 orang memilih Soda?
Perhitungan Manual di R
# Jumlah total individu
n <- 15
# Frekuensi tiap kategori
x <- c(6, 4, 3, 2) # T, K, J, S
# Probabilitas preferensi
p <- c(0.30, 0.25, 0.25, 0.20)
# Hitung koefisien multinomial
faktorial_total <- factorial(n)
faktorial_x <- prod(factorial(x))
koefisien <- faktorial_total / faktorial_x
# Hitung probabilitas gabungan
peluang <- koefisien * prod(p^x)
peluang
## [1] 0.01122386
Interpretasi: Probabilitas bahwa secara acak akan diperoleh 6 orang memilih Teh, 4 orang Kopi, 3 orang Jus, dan 2 orang Soda (dengan asumsi probabilitas preferensi masing-masing adalah 0.30, 0.25, 0.25, dan 0.20) adalah sebesar 0.01122386. Nilai ini dihitung menggunakan rumus distribusi multinomial, yang merupakan generalisasi dari distribusi binomial ketika terdapat lebih dari dua kategori yang bersifat saling eksklusif.
Model ini digunakan untuk menjelaskan hubungan antara satu variabel respons kategorik (memiliki lebih dari dua kategori) dan satu atau lebih variabel prediktor.
Misalkan \(Y\) terdiri atas \(K\) kategori, dan kategori referensi (baseline) adalah kategori \(K\). Maka bentuk umum model logit multinomial:
\[ \log\left(\frac{P(Y = j)}{P(Y = K)}\right) = \alpha_j + \beta_{j1}x_1 + \beta_{j2}x_2 + \dots + \beta_{jp}x_p \]
Untuk \(j = 1, 2, \dots, K-1\)
Model baseline-category logit merupakan variasi dari regresi logistik yang diterapkan pada variabel respons dengan lebih dari dua kategori. Dalam pendekatan ini, satu kategori digunakan sebagai acuan (baseline), kemudian kategori lainnya dibandingkan terhadap kategori baseline tersebut.
Rumus modelnya:
\[ \log\left(\frac{\pi_j}{\pi_K}\right) = \eta_j, \quad j = 1,\dots,K-1 \]
dengan:
Model ini akan menghasilkan \(K - 1\) persamaan logit.
Catatan: Secara default, multinom()
dari package nnet
di R menjadikan kategori terakhir sebagai
baseline. Namun, kita juga bisa menentukan kategori referensi secara
manual.
Model Logit dengan Satu Prediktor
Jika hanya ada satu variabel prediktor, misalkan \(x\), maka bentuk modelnya menjadi:
\[ \log\left(\frac{\pi_j}{\pi_K}\right) = \alpha_j + \beta_j x, \quad j = 1,\dots,K-1 \] Contoh Kasus: Respons 3 Kategori
Misalkan variabel respons \(Y\) memiliki 3 kategori: \(Y \in \{A, B, C\}\), dan kita jadikan kategori \(C\) sebagai baseline. Maka modelnya:
\[ \log\left(\frac{\pi_A}{\pi_C}\right) = \alpha_1 + \beta_1 x \\ \log\left(\frac{\pi_B}{\pi_C}\right) = \alpha_2 + \beta_2 x \] Relasi Antar Kategori Jika kita ingin mengetahui logit antara dua kategori non-baseline, misalkan antara kategori 1 (\(\pi_1\)) dan kategori 2 (\(\pi_2\)), maka kita bisa menggunakan rumus berikut:
\[ \log\left(\frac{\pi_1}{\pi_2}\right) = \log\left(\frac{\pi_1 / \pi_3}{\pi_2 / \pi_3}\right) = \log\left(\frac{\pi_1}{\pi_3}\right) - \log\left(\frac{\pi_2}{\pi_3}\right) \]
Jika model logit tiap kategori terhadap baseline (\(\pi_3\)) adalah:
\[ \log\left(\frac{\pi_1}{\pi_3}\right) = \alpha_1 + \beta_1 x \\ \log\left(\frac{\pi_2}{\pi_3}\right) = \alpha_2 + \beta_2 x \]
Maka logit antara kategori 1 dan 2 menjadi:
\[ \log\left(\frac{\pi_1}{\pi_2}\right) = (\alpha_1 - \alpha_2) + (\beta_1 - \beta_2)x \] Ringkasan Model Logit Baseline-Kategori
Catatan Implementasi di R:
Fungsi multinom()
dari paket nnet
digunakan
untuk membangun model ini. Pengaturan kategori baseline dapat diubah
menggunakan fungsi relevel()
.
Untuk menentukan parameter model, digunakan metode Maximum Likelihood Estimation (MLE). Estimasi dilakukan melalui algoritma iteratif seperti Newton-Raphson.
Fungsi log-likelihood yang dimaksimalkan:
\[ \ell(\boldsymbol{\beta}) = \sum_{i=1}^{n} \sum_{j=1}^{K} y_{ij} \log(\pi_{ij}) \]
di mana:
Sebuah perusahaan retail nasional ingin mengevaluasi preferensi staf gudang terhadap jenis kendaraan operasional yang digunakan untuk distribusi barang, yaitu antara: Sepeda Motor, Mobil Box, dan Pick-up.
Untuk itu, dilakukan survei terhadap 180 karyawan dengan mencatat beberapa informasi sebagai berikut:
Tujuan Penelitian
Mengetahui apakah unit kerja, usia, dan pengalaman kerja memiliki keterkaitan dalam preferensi pemilihan kendaraan operasional.
# Penyiapan pustaka
library(nnet)
## Warning: package 'nnet' was built under R version 4.3.3
set.seed(999)
# Jumlah responden
JUMLAH_RESP <- 180
# Simulasi variabel DIVISI
DIVISI <- sample(c("Logistik", "Distribusi", "Stok Barang"), size = JUMLAH_RESP, replace = TRUE)
# Simulasi usia dengan distribusi normal
USIA_KERJA <- round(rnorm(JUMLAH_RESP, mean = 37, sd = 6))
# Simulasi pengalaman kerja (tahun)
LAMA_KERJA <- round(pmax(rnorm(JUMLAH_RESP, mean = 5, sd = 2.5), 0)) # hanya nilai positif
# Simulasi preferensi kendaraan berdasarkan DIVISI
KENDARAAN <- sapply(DIVISI, function(bidang) {
if (bidang == "Logistik") {
sample(c("Motorcycle", "Van", "Pickup"), size = 1, prob = c(0.3, 0.5, 0.2))
} else if (bidang == "Distribusi") {
sample(c("Motorcycle", "Van", "Pickup"), size = 1, prob = c(0.2, 0.3, 0.5))
} else {
sample(c("Motorcycle", "Van", "Pickup"), size = 1, prob = c(0.4, 0.4, 0.2))
}
})
# Gabungkan semua ke dalam satu data.frame
DATA_RESP <- data.frame(
VEHICLE = factor(KENDARAAN),
AGE = USIA_KERJA,
DIVISION = factor(DIVISI),
SENIORITY = LAMA_KERJA
)
# Atur baseline kategori untuk VEHICLE
DATA_RESP$VEHICLE <- relevel(DATA_RESP$VEHICLE, ref = "Pickup")
# Tampilkan 10 data pertama
head(DATA_RESP, 10)
## VEHICLE AGE DIVISION SENIORITY
## 1 Van 34 Stok Barang 8
## 2 Van 38 Logistik 4
## 3 Van 24 Stok Barang 8
## 4 Van 26 Logistik 8
## 5 Van 44 Distribusi 5
## 6 Pickup 38 Logistik 3
## 7 Van 39 Distribusi 7
## 8 Motorcycle 38 Distribusi 9
## 9 Pickup 43 Distribusi 6
## 10 Pickup 38 Stok Barang 6
#Visualisasi
library(ggplot2)
ggplot(DATA_RESP, aes(x = VEHICLE, fill = DIVISION)) +
geom_bar(position = "dodge") +
labs(
title = "Distribusi Pemilihan Kendaraan Operasional",
x = "Jenis Kendaraan", y = "Jumlah Pegawai"
) +
theme_minimal()
Setelah data disimulasikan, langkah berikutnya adalah melakukan
estimasi model menggunakan pendekatan baseline-category
logit. Model ini mengevaluasi pengaruh beberapa prediktor
terhadap kategori pilihan kendaraan, dengan
Pickup
sebagai baseline.
# Estimasi model logit multinomial
model_kendaraan <- multinom(VEHICLE ~ AGE + DIVISION + SENIORITY, data = DATA_RESP)
## # weights: 18 (10 variable)
## initial value 197.750212
## iter 10 value 184.287340
## final value 184.162668
## converged
# Ringkasan hasil model
summary(model_kendaraan)
## Call:
## multinom(formula = VEHICLE ~ AGE + DIVISION + SENIORITY, data = DATA_RESP)
##
## Coefficients:
## (Intercept) AGE DIVISIONLogistik DIVISIONStok Barang
## Motorcycle 0.8709305 -0.04393051 1.503878 1.7930247
## Van 0.9033808 -0.02814743 1.260367 0.6306604
## SENIORITY
## Motorcycle -0.08651643
## Van -0.02832544
##
## Std. Errors:
## (Intercept) AGE DIVISIONLogistik DIVISIONStok Barang
## Motorcycle 1.444861 0.03631097 0.5307996 0.5123638
## Van 1.292663 0.03227399 0.4384735 0.4590667
## SENIORITY
## Motorcycle 0.08579943
## Van 0.07757530
##
## Residual Deviance: 368.3253
## AIC: 388.3253
Untuk mengevaluasi apakah masing-masing prediktor berkontribusi signifikan secara statistik, kita hitung nilai z-score dan p-value berdasarkan rasio antara koefisien dan standard error.
# Hitung nilai Z-statistik
Z_stat <- summary(model_kendaraan)$coefficients / summary(model_kendaraan)$standard.errors
# Hitung p-value 2 sisi
p_value <- 2 * (1 - pnorm(abs(Z_stat)))
# Bulatkan hasil
round(p_value, 4)
## (Intercept) AGE DIVISIONLogistik DIVISIONStok Barang SENIORITY
## Motorcycle 0.5467 0.2263 0.0046 0.0005 0.3133
## Van 0.4846 0.3831 0.0040 0.1695 0.7150
Berdasarkan hasil perhitungan z-score dan p-value, dapat disimpulkan bahwa beberapa prediktor memiliki pengaruh yang signifikan terhadap preferensi kendaraan operasional pegawai.
Setelah model multinomial berhasil diestimasi, langkah selanjutnya adalah melakukan prediksi kategori kendaraan untuk masing-masing observasi. Prediksi ini kemudian divalidasi dengan membandingkan hasil klasifikasi dengan kategori aktual dari data simulasi.
# Prediksi kategori kendaraan untuk seluruh observasi
DATA_RESP$HASIL_PREDIKSI <- predict(model_kendaraan)
# Buat tabel perbandingan antara aktual dan prediksi
table(TERPREDIKSI = DATA_RESP$HASIL_PREDIKSI, SEBENARNYA = DATA_RESP$VEHICLE)
## SEBENARNYA
## TERPREDIKSI Pickup Motorcycle Van
## Pickup 28 9 22
## Motorcycle 11 23 18
## Van 17 18 34
Berdasarkan hasil estimasi dan validasi model regresi logistik multinomial terhadap data simulasi preferensi kendaraan operasional, ditemukan bahwa divisi kerja pegawai merupakan faktor yang paling berpengaruh dalam menentukan pilihan kendaraan (Motorcycle atau Van) dibandingkan dengan referensi Pickup. Hasil ini diperkuat oleh signifikansi statistik dari variabel DIVISION pada kedua kategori tersebut, sementara usia dan pengalaman kerja tidak menunjukkan pengaruh yang signifikan. Melalui model ini, perusahaan dapat memetakan kecenderungan pemilihan kendaraan berdasarkan latar belakang unit kerja, yang dapat dimanfaatkan untuk merencanakan distribusi logistik atau pengadaan kendaraan operasional yang lebih tepat sasaran sesuai dengan kebutuhan karyawan di masing-masing divisi.
Ketika variabel respons \(Y\) memiliki sifat urutan (ordinal), seperti: tingkat kesadaran kesehatan (Rendah, Sedang, Tinggi), maka regresi logistik ordinal menjadi pendekatan yang tepat untuk menganalisis hubungan antara \(Y\) dan satu atau lebih prediktor.
Karakteristik:
Salah satu pendekatan paling umum dalam regresi ordinal adalah Cumulative Logit Model, dengan asumsi dasar proportional odds (peluang proporsional).
Bentuk model:
\[ \log\left( \frac{P(Y \le j)}{P(Y > j)} \right) = \alpha_j + \boldsymbol{\beta} x \] Penjelasan:
Jika terdapat \(K\) kategori, maka akan terbentuk sebanyak \((K - 1)\) model logit kumulatif.
Dalam regresi logistik ordinal (cumulative logit model), koefisien regresi (\(\beta\)) menunjukkan arah dan besar pengaruh prediktor terhadap probabilitas kumulatif dari suatu kategori dan yang lebih rendah darinya.
Prinsip Umum Interpretasi Koefisien
Model membandingkan:
\[ \log\left(\frac{P(Y \le j)}{P(Y >
j)}\right) = \alpha_j + \beta x \]
untuk setiap batas kumulatif antar kategori.
Jika \(\beta > 0\): peningkatan nilai \(x\) akan meningkatkan peluang berada di kategori bawah atau menengah (dan menurunkan peluang kategori tertinggi).
Jika \(\beta < 0\): peningkatan nilai \(x\) justru menaikkan kemungkinan berada di kategori lebih tinggi.
Odds Ratio (Rasio Peluang)
Koefisien juga dapat diterjemahkan ke dalam bentuk Odds Ratio (OR):
\[ \text{OR} = e^{\beta} \]
Dalam studi ini, sebuah survei dilakukan terhadap 250 warga di beberapa kota besar. Mereka diminta menilai kepuasan terhadap kebersihan taman umum dengan skala ordinal sebagai berikut:
1 = Tidak Bersih
2 = Lumayan
3 = Sangat Bersih
Salah satu prediktor yang dianalisis adalah jumlah petugas kebersihan yang aktif di taman tersebut per hari.
set.seed(456)
JUMLAH_RESPONDEN <- 250
PETUGAS <- sample(1:10, JUMLAH_RESPONDEN, replace = TRUE)
KONDISI <- cut(
1.2*PETUGAS + rnorm(JUMLAH_RESPONDEN, sd = 2),
breaks = c(-Inf, 6, 10, Inf),
labels = c("Tidak Bersih", "Lumayan", "Sangat Bersih"),
ordered_result = TRUE
)
DATA_TAMAN <- data.frame(KONDISI, PETUGAS)
head(DATA_TAMAN)
## KONDISI PETUGAS
## 1 Tidak Bersih 5
## 2 Tidak Bersih 3
## 3 Lumayan 6
## 4 Lumayan 5
## 5 Lumayan 9
## 6 Sangat Bersih 10
Kita akan membentuk model regresi logistik ordinal untuk mengetahui pengaruh jumlah petugas kebersihan terhadap persepsi warga.
library(MASS)
MODEL_KOTA <- polr(KONDISI ~ PETUGAS, data = DATA_TAMAN, Hess = TRUE)
summary(MODEL_KOTA)
## Call:
## polr(formula = KONDISI ~ PETUGAS, data = DATA_TAMAN, Hess = TRUE)
##
## Coefficients:
## Value Std. Error t value
## PETUGAS 1.087 0.1039 10.47
##
## Intercepts:
## Value Std. Error t value
## Tidak Bersih|Lumayan 5.5074 0.5864 9.3911
## Lumayan|Sangat Bersih 8.8198 0.8460 10.4252
##
## Residual Deviance: 272.8232
## AIC: 278.8232
Langkah ini bertujuan untuk menilai signifikansi statistik dari prediktor dan intercept. Kita gunakan distribusi normal standar karena model polr() tidak langsung menyediakan nilai p-value.
## ----p_value_calc--------------------------------------------------------
TABEL_KOEF <- coef(summary(MODEL_KOTA))
## Hitung p-value dari t-value
PVAL <- 2 * pnorm(-abs(TABEL_KOEF[, "t value"]))
TABEL_KOEF <- cbind(TABEL_KOEF, "p-value" = round(PVAL, 4))
TABEL_KOEF
## Value Std. Error t value p-value
## PETUGAS 1.087161 0.1038666 10.466894 0
## Tidak Bersih|Lumayan 5.507367 0.5864440 9.391121 0
## Lumayan|Sangat Bersih 8.819769 0.8460009 10.425247 0
## ----predicted_probabilities--------------------------------------------
DATA_BARU <- data.frame(PETUGAS = 4:8)
prediksi_prob <- predict(MODEL_KOTA, newdata = DATA_BARU, type = "probs")
round(prediksi_prob, 4)
## Tidak Bersih Lumayan Sangat Bersih
## 1 0.7611 0.2276 0.0113
## 2 0.5179 0.4493 0.0328
## 3 0.2659 0.6427 0.0914
## 4 0.1088 0.6614 0.2298
## 5 0.0395 0.4910 0.4694
## ----plot_probabilities-------------------------------------------------
library(tidyr)
## Warning: package 'tidyr' was built under R version 4.3.3
library(ggplot2)
PROB_TBL <- as.data.frame(prediksi_prob)
PROB_TBL$PETUGAS <- 4:8
PROB_TBL_LONG <- pivot_longer(PROB_TBL, cols = -PETUGAS, names_to = "KATEGORI", values_to = "PROB")
ggplot(PROB_TBL_LONG, aes(x = PETUGAS, y = PROB, color = KATEGORI)) +
geom_line(size = 1.2) +
labs(title = "Perubahan Probabilitas terhadap Jumlah Petugas",
x = "Jumlah Petugas Kebersihan",
y = "Probabilitas",
color = "Tingkat Kebersihan") +
theme_minimal()
Interpretasi
Grafik menunjukkan perubahan probabilitas tiga kategori tingkat kebersihan taman kota berdasarkan jumlah petugas kebersihan.
Ketika jumlah petugas rendah (4–5 orang), probabilitas taman dinilai Tidak Bersih (garis biru) sangat tinggi, sementara probabilitas Sangat Bersih (garis hijau) sangat rendah.
Saat jumlah petugas meningkat, probabilitas persepsi “Tidak Bersih” menurun tajam, sedangkan kategori Lumayan (merah) naik, dan kemudian kategori “Sangat Bersih” meningkat signifikan setelah petugas melewati 6 orang.
Pada jumlah petugas ≥8, Sangat Bersih menjadi persepsi yang paling dominan.
Artinya: Peningkatan jumlah petugas kebersihan secara nyata meningkatkan peluang masyarakat menilai taman menjadi lebih bersih.
Model cumulative logit (seperti yang digunakan dalam polr) mengasumsikan bahwa efek dari prediktor tetap konstan di seluruh batas antar kategori. Artinya, perubahan nilai prediktor (jumlah petugas) dianggap memberi pengaruh yang sama antara:
Jika asumsi ini tidak dipenuhi, maka hasil interpretasi bisa menjadi bias. Salah satu solusi adalah mempertimbangkan model yang tidak mengasumsikan proportional odds.
Jika asumsi cumulative logit tidak sesuai, kita bisa pertimbangkan: - Adjacent-category logit: membandingkan kategori yang berurutan (misalnya: “Tidak Bersih” vs “Lumayan”, lalu “Lumayan” vs “Sangat Bersih”). - Continuation-ratio logit: digunakan ketika data memiliki sifat bertahap dan bersifat urutan seperti progresi penyakit. - Generalized Ordered Logit: fleksibel karena membolehkan koefisien berbeda di setiap cut-off.
Berdasarkan hasil analisis regresi logistik ordinal, diperoleh bahwa jumlah petugas kebersihan (PETUGAS) memiliki koefisien positif sebesar 1.087 dengan nilai t sangat tinggi (10.47) dan nilai p-value ≈ 0, menunjukkan bahwa prediktor ini sangat signifikan secara statistik dalam memengaruhi persepsi warga terhadap kebersihan taman. Intercept antar kategori, yaitu transisi dari “Tidak Bersih” ke “Lumayan” (5.51) dan “Lumayan” ke “Sangat Bersih” (8.82), juga signifikan. Hal ini mengindikasikan bahwa peningkatan satu petugas kebersihan secara substansial menaikkan peluang taman dinilai lebih bersih oleh warga. Secara keseluruhan, model mendukung hipotesis bahwa jumlah petugas berpengaruh positif terhadap persepsi kebersihan taman kota, dan dapat digunakan untuk perencanaan distribusi tenaga kebersihan secara lebih efisien.
Model regresi logistik ordinal yang umum dipakai adalah model Logit Kumulatif (Cumulative Logit Model) dengan dasar asumsi Odds Proporsional.
Apa Itu Asumsi Paralel? Asumsi paralel (juga dikenal sebagai parallel lines assumption) mengharuskan bahwa pengaruh prediktor (β) harus identik di seluruh batas kategori pada variabel respon ordinal.
Bentuk model matematis umum:
\[ \log\left( \frac{P(Y \leq j)}{P(Y > j)} \right) = \alpha_j + \beta x \]
Untuk \(j = 1, \dots, c-1\), dengan ketentuan:
\(\alpha_j\): nilai konstanta (intercept) yang berbeda untuk tiap ambang batas.
\(\beta\): koefisien yang tetap sama di semua persamaan.
Ilustrasi Visual
Jika asumsi ini benar, garis logit kumulatif dari tiap kategori terhadap variabel prediktor akan paralel dan hanya berbeda dari sisi posisi awal (intersep), bukan dari kemiringan garisnya.
Apa yang Terjadi Jika Asumsi Tidak Dipenuhi?
Jika asumsi ini tidak berlaku:
Efek variabel bebas bisa berubah antar kategori.
Model menjadi bias dan tidak sah.
Harus menggunakan pendekatan lain seperti:
Generalized Ordinal Logistic Model
Partial Proportional Odds Model
Pengujian Asumsi Paralel
Untuk mengevaluasi apakah asumsi ini terpenuhi, dapat dilakukan: -
Uji Rasio Likelihood (LRT) antara model proporsional
dan non-proporsional. - Uji Brant, yang tersedia di
paket brant
pada R
Contoh Implematansi di R
library(MASS)
library(brant)
## Warning: package 'brant' was built under R version 4.3.3
# Simulasi data dan fitting model ordinal
set.seed(321)
dat <- data.frame(
y_ord = factor(sample(c("Rendah", "Sedang", "Tinggi"), 150, replace = TRUE), ordered = TRUE),
x1 = rnorm(150),
x2 = runif(150)
)
# Model ordinal logit
model_ordinal <- polr(y_ord ~ x1 + x2, data = dat, Hess = TRUE)
# Uji Brant
brant(model_ordinal)
## --------------------------------------------
## Test for X2 df probability
## --------------------------------------------
## Omnibus 1.67 2 0.43
## x1 0.01 1 0.91
## x2 1.65 1 0.2
## --------------------------------------------
##
## H0: Parallel Regression Assumption holds
Ringkasan
Asumsi paralel memainkan peranan penting dalam memastikan validitas model logit kumulatif. Asumsi ini menyederhanakan model dengan menganggap efek prediktor konsisten pada tiap ambang kategori. Jika ternyata tidak terpenuhi, maka kita harus mempertimbangkan pendekatan model yang lebih fleksibel agar interpretasi dan hasil analisis tetap akurat dan dapat diandalkan.
Dalam bidang statistika terapan, analisis data kategorik memainkan peran penting dalam memahami fenomena-fenomena yang umum terjadi di berbagai bidang seperti sosiologi, kesehatan masyarakat, pemasaran, dan pendidikan. Misalnya, data tentang jenis kelamin, kategori usia, status pekerjaan, atau preferensi pelanggan biasanya disimpan dalam bentuk kategori, bukan angka kontinu.
Beberapa metode yang umum digunakan untuk menganalisis data jenis ini antara lain:
Tabel Kontingensi
Model Log-Linear
Model Regresi Logistik
Setiap metode memiliki fungsi dan kelebihan masing-masing tergantung pada tujuan analisis yang ingin dicapai.
Tabel Kontingensi
Tabel kontingensi merupakan langkah awal yang sering digunakan untuk mengeksplorasi hubungan antara dua atau lebih variabel diskrit. Misalnya, kita ingin mengetahui hubungan antara paparan obat tertentu dan hasil pengobatan (sembuh/tidak). Tabel ini membantu mengidentifikasi pola distribusi data dan potensi asosiasi dengan menggunakan uji seperti chi-square.
Model Log-Linear
Berbeda dengan regresi logistik, model log-linear tidak membedakan secara eksplisit antara variabel dependen dan independen. Sebaliknya, pendekatan ini memperlakukan semua variabel kategori secara simetris dan cocok digunakan saat tujuan analisis adalah menggambarkan struktur interaksi antar variabel. Model ini termasuk dalam keluarga Generalized Linear Models (GLM) dengan asumsi bahwa data mengikuti distribusi Poisson.
Struktur model log-linear dirancang untuk mengukur efek utama dan efek interaksi antar variabel. Sebagai contoh, dalam sebuah studi kesehatan yang melibatkan variabel jenis kelamin, status merokok, dan kejadian kanker paru-paru, model log-linear bisa menunjukkan apakah ada interaksi yang signifikan antara variabel-variabel tersebut dalam memengaruhi kejadian kanker. Uji Likelihood Ratio Test (LRT) sering digunakan untuk membandingkan model dengan dan tanpa interaksi.
Regresi Logistik
Sementara model log-linear lebih fokus pada struktur asosiasi, regresi logistik bertujuan untuk prediksi. Model ini cocok saat terdapat satu variabel respon biner (misalnya: terkena penyakit, ya/tidak), dan satu atau lebih variabel prediktor baik berupa kategori maupun numerik. Dengan menggunakan fungsi logit, model ini memodelkan probabilitas kejadian outcome tertentu berdasarkan nilai dari prediktor.
Regresi logistik sangat umum digunakan dalam studi epidemiologi dan ilmu sosial, dan dapat dikembangkan menjadi versi multinomial (untuk lebih dari dua kategori) maupun ordinal (untuk kategori berurutan).
Ringkasan Pendekatan
Secara umum, pendekatan-pendekatan yang sering digunakan dalam analisis data kategorik antara lain:
Tabel Kontingensi
Digunakan untuk mendeskripsikan distribusi frekuensi dua atau lebih
variabel kategorik secara tabular.
Model Log-Linear
Digunakan untuk mengeksplorasi pola interaksi antara variabel-variabel
kategori tanpa menetapkan satu variabel sebagai target.
Regresi Logistik
Digunakan ketika fokusnya adalah memodelkan probabilitas suatu outcome
berdasarkan prediktor, dengan membedakan antara variabel respon dan
prediktor.
Meskipun ketiganya dapat digunakan pada data kategorik, pendekatan dan interpretasinya sangat berbeda.
Tabel Kontingensi
Tabel kontingensi menyajikan distribusi frekuensi dari gabungan dua atau lebih kategori, memungkinkan kita mengevaluasi hubungan antara dua variabel nominal atau ordinal.
Contoh tabel 2x2:
table_kat <- matrix(c(25, 35, 45, 55), nrow = 2, byrow = TRUE,
dimnames = list(
Terapi = c("Aspirin", "Kontrol"),
Status = c("Sembuh", "Tidak")
))
table_kat
## Status
## Terapi Sembuh Tidak
## Aspirin 25 35
## Kontrol 45 55
Model Loglinear
Model log-linear digunakan ketika tujuan kita adalah mengidentifikasi struktur asosiasi atau ketergantungan antar kategori dalam tabel, tanpa membedakan antara respon dan prediktor.
Model ini memodelkan log dari ekspektasi frekuensi tiap sel:
\[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \] Model ini cocok ketika: - Kita ingin menyelidiki hubungan antar kategori. - Tidak ada variabel respon eksplisit. - Digunakan dalam tabel multi-arah (misal 2x2x2).
library(MASS)
loglm(~ Terapi + Status, data = table_kat)
## Call:
## loglm(formula = ~Terapi + Status, data = table_kat)
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0.1695917 1 0.6804750
## Pearson 0.1693122 1 0.6807238
Model Regresi Logistik Model regresi logistik digunakan ketika kita memiliki satu variabel target yang bersifat kategorik, khususnya biner.
Bentuk model: \[ \log\left( \frac{p}{1 - p} \right) = \beta_0 + \beta_1 x \]
Contoh : Misalnya kita ingin memodelkan kemungkinan seseorang sembuh berdasarkan apakah mereka menerima terapi:
# Data simulasi
df <- data.frame(
Sembuh = c(1, 1, 0, 0, 1, 0, 1, 0),
Terapi = c(1, 1, 1, 0, 0, 0, 1, 0)
)
# Regresi logistik
model_logistik <- glm(Sembuh ~ Terapi, data = df, family = binomial)
summary(model_logistik)
##
## Call:
## glm(formula = Sembuh ~ Terapi, family = binomial, data = df)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.099 1.155 -0.951 0.341
## Terapi 2.197 1.633 1.346 0.178
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 11.0904 on 7 degrees of freedom
## Residual deviance: 8.9974 on 6 degrees of freedom
## AIC: 12.997
##
## Number of Fisher Scoring iterations: 4
Perbandingan dari Ketiga Pendekatan
Aspek | Tabel Kontingensi | Model Loglinear | Regresi Logistik |
---|---|---|---|
Tujuan | Deskripsi frekuensi | Deteksi asosiasi | Prediksi probabilitas |
Variabel Dependen | Tidak ada | Tidak ada (simetris) | Ada (eksplisit) |
Distribusi | Tidak diasumsikan | Poisson (frekuensi sel) | Binomial (probabilitas) |
Bentuk Model | Tidak ada | GLM: log(μ) ~ efek | GLM: logit(p) ~ prediktor |
Cocok untuk | Eksplorasi awal | Tabel > 2 variabel | Studi prediktif |
Tabel kontingensi digunakan untuk menyajikan data frekuensi dari berbagai kombinasi kategori dua atau lebih variabel. Sebagai contoh, kita ingin melihat hubungan antara jenis perlakuan dan hasil kejadian pada pasien.
Berikut contoh sederhana tabel 2x2:
# Tabel kontingensi 2x2
tabel_kategori <- matrix(c(28, 52, 23, 57), nrow = 2,
dimnames = list(
Perlakuan = c("Ekstrak_A", "Kontrol"),
Respon = c("Positif", "Negatif")
))
tabel_kategori
## Respon
## Perlakuan Positif Negatif
## Ekstrak_A 28 23
## Kontrol 52 57
Model log-linear pada tabel kontingensi memodelkan nilai ekspektasi dari frekuensi sel berdasarkan efek utama dan interaksi:
Model log-linier untuk tabel I x J dapat dituliskan: \[ \log(\mu_{ij}) = \mu + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \]
Model ini sangat sesuai untuk menyelidiki hubungan antar variabel kategori tanpa mendefinisikan mana yang sebagai prediktor dan mana respon.
Model saturated (atau model penuh) adalah model log-linear yang menyertakan semua efek utama dan interaksi tingkat dua, sehingga cocok sempurna terhadap data.
Karakteristik:
Tidak mengasumsikan adanya independensi antar variabel.
Memberikan kecocokan sempurna terhadap frekuensi yang diamati.
Simulasi Tabel dan Model Berikut ini contoh implementasi model saturated menggunakan data 2x2:
library(MASS)
# Data simulasi baru
data_saturasi <- matrix(c(33, 47, 36, 54), nrow = 2, byrow = TRUE,
dimnames = list(
Obat = c("ObatX", "Plasebo"),
Status = c("Ya", "Tidak")
))
# Tabel sebagai input
ftable(data_saturasi)
## Status Ya Tidak
## Obat
## ObatX 33 47
## Plasebo 36 54
Model saturated dapat dibangun dengan fungsi loglm() dari package MASS, dengan menyertakan interaksi penuh antar variabel:
# Model saturated (dengan interaksi penuh)
model_saturasi <- loglm(~ Obat * Status, data = data_saturasi)
# Ringkasan hasil
summary(model_saturasi)
## Formula:
## ~Obat * Status
## attr(,"variables")
## list(Obat, Status)
## attr(,"factors")
## Obat Status Obat:Status
## Obat 1 0 1
## Status 0 1 1
## attr(,"term.labels")
## [1] "Obat" "Status" "Obat:Status"
## attr(,"order")
## [1] 1 1 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0 0 1
## Pearson 0 0 1
Model saturated akan selalu memberikan deviance = 0 dan p-value = 1, karena tidak menyisakan variansi yang tidak dijelaskan (fit sempurna).
Model independen mengasumsikan tidak ada interaksi antara dua variabel kategorik. Artinya, kombinasi antar kategori tidak memberikan informasi lebih terhadap satu sama lain.
Modelnya dituliskan sebagai: \[ \log(\mu_{ij}) = \mu + \lambda_i^A + \lambda_j^B \] Tujuan Model:
library(MASS)
# Simulasi data 2x2
data_indep <- matrix(c(29, 41, 38, 52), nrow = 2, byrow = TRUE,
dimnames = list(
Perlakuan = c("Terapi_Z", "Kontrol"),
Respon = c("Ya", "Tidak")
))
# Model independen (tanpa interaksi)
model_indep <- loglm(~ Perlakuan + Respon, data = data_indep)
# Ringkasan hasil
summary(model_indep)
## Formula:
## ~Perlakuan + Respon
## attr(,"variables")
## list(Perlakuan, Respon)
## attr(,"factors")
## Perlakuan Respon
## Perlakuan 1 0
## Respon 0 1
## attr(,"term.labels")
## [1] "Perlakuan" "Respon"
## attr(,"order")
## [1] 1 1
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0.01019208 1 0.9195855
## Pearson 0.01018971 1 0.9195949
Interpretasi
Nilai statistik Chi-Square dan p-value akan menunjukkan apakah asumsi independensi layak dipertahankan.
Jika p-value < 0.05, maka interaksi perlu dipertimbangkan (data tidak independen).
Odds ratio (OR) adalah ukuran asosiasi yang umum digunakan dalam tabel kontingensi 2x2. OR mengukur peluang relatif dari suatu kejadian antara dua kelompok.
\[ OR = \frac{n_{11} n_{22}}{n_{12} n_{21}} \] Dimana :
\(n_{11}\) : Frekuensi pada sel baris pertama kolom pertama
\(n_{12}\) : Frekuensi pada sel baris pertama kolom kedua
\(n_{21}\) : Frekuensi pada sel baris kedua kolom pertama
\(n_{22}\) : Frekuensi pada sel baris kedua kolom kedua
Interprestasi Nilai OR:
\(OR = 1\) : Tidak ada hubungan (asosiasi netral)
\(OR > 1\) : Hubungan positif (kecenderungan kejadian lebih tinggi pada grup pertama)
\(OR < 1\) : Hubungan negatif (kecenderungan kejadian lebih rendah pada grup pertama)
Contoh Perhitungan OR
# Data 2x2
data_or <- matrix(c(35, 40, 30, 55), nrow = 2, byrow = TRUE)
# Hitung OR manual
OR_manual <- (data_or[1,1] * data_or[2,2]) / (data_or[1,2] * data_or[2,1])
OR_manual
## [1] 1.604167
Dalam model log-linear, khususnya model saturated, estimasi parameter bertujuan untuk mengkuantifikasi efek utama dan interaksi antar variabel kategori. Estimasi dilakukan berdasarkan asumsi:
Salah satu ukuran yang bisa diekstrak dari tabel kontingensi adalah Odds Ratio (OR) dan log-odds, yang secara langsung berkaitan dengan parameter dalam model log-linear 2x2.
Rumus estimasi log-odds:
\[ \log(OR) = \log\left( \frac{n_{11} \cdot n_{22}}{n_{12} \cdot n_{21}} \right) \] Contoh perhitungan:
# Estimasi log(OR)
logOR <- log((data_or[1,1] * data_or[2,2]) / (data_or[1,2] * data_or[2,1]))
logOR
## [1] 0.4726044
Untuk menilai apakah model yang lebih sederhana (tanpa interaksi) cukup baik, kita bisa melakukan Likelihood Ratio Test (LRT) atau uji deviasi (G²).
Tujuan Pengujuan:
Menentukan apakah interaksi antar variabel diperlukan.
Jika model sederhana cukup baik (p-value > 0.05), maka model saturated dianggap tidak memberikan informasi tambahan yang signifikan.
# Bandingkan model independen vs saturated
anova(model_indep, model_saturasi)
## LR tests for hierarchical log-linear models
##
## Model 1:
## ~Perlakuan + Respon
## Model 2:
## ~Obat * Status
##
## Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1 0.01019208 1
## Model 2 0.00000000 0 0.01019208 1 0.91959
## Saturated 0.00000000 0 0.00000000 0 1.00000
Interpretasi
Nilai deviance selisih (Delta(Dev)) = 0.01019, artinya hanya ada sedikit perbedaan kecocokan model antara model independen dan saturated.
P-value = 0.91959 sangat tinggi (≫ 0.05), sehingga tidak ada bukti cukup untuk menolak hipotesis nol (bahwa model sederhana sudah cukup).
Artinya, penambahan interaksi dalam model saturated tidak memberikan peningkatan signifikan dalam kecocokan model terhadap data.
Kesimpulan: Kita tidak perlu menggunakan model saturated, cukup memakai model independen (lebih sederhana dan interpretatif).
Data ini berasal dari survei General Social Survey (GSS) tahun 1989, yang mengklasifikasikan responden berdasarkan pandangan politik (liberal, moderat, konservatif) dan pilihan mereka dalam pemilihan presiden 1988 (Dukakis atau Bush).:contentReferenceoaicite:6
Sumber: University of Illinois - Log-linear Models for Contingency Tables
# Membuat tabel kontingensi
tabel_politik <- matrix(c(197, 65, 148, 186, 68, 242), nrow = 3, byrow = TRUE,
dimnames = list(
Pandangan = c("Liberal", "Moderat", "Konservatif"),
Pilihan = c("Dukakis", "Bush")
))
tabel_politik
## Pilihan
## Pandangan Dukakis Bush
## Liberal 197 65
## Moderat 148 186
## Konservatif 68 242
# Model saturated
model_saturated <- loglm(~ Pandangan * Pilihan, data = tabel_politik)
# Ringkasan model
summary(model_saturated)
## Formula:
## ~Pandangan * Pilihan
## attr(,"variables")
## list(Pandangan, Pilihan)
## attr(,"factors")
## Pandangan Pilihan Pandangan:Pilihan
## Pandangan 1 0 1
## Pilihan 0 1 1
## attr(,"term.labels")
## [1] "Pandangan" "Pilihan" "Pandangan:Pilihan"
## attr(,"order")
## [1] 1 1 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0 0 1
## Pearson 0 0 1
# Model independen
model_independen <- loglm(~ Pandangan + Pilihan, data = tabel_politik)
# Ringkasan model
summary(model_independen)
## Formula:
## ~Pandangan + Pilihan
## attr(,"variables")
## list(Pandangan, Pilihan)
## attr(,"factors")
## Pandangan Pilihan
## Pandangan 1 0
## Pilihan 0 1
## attr(,"term.labels")
## [1] "Pandangan" "Pilihan"
## attr(,"order")
## [1] 1 1
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 170.4841 2 0
## Pearson 162.6962 2 0
# Perbandingan model
anova(model_independen, model_saturated)
## LR tests for hierarchical log-linear models
##
## Model 1:
## ~Pandangan + Pilihan
## Model 2:
## ~Pandangan * Pilihan
##
## Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1 170.4841 2
## Model 2 0.0000 0 170.4841 2 0
## Saturated 0.0000 0 0.0000 0 1
# Menghitung odds ratio antara Liberal dan Konservatif
OR <- (tabel_politik["Liberal", "Dukakis"] * tabel_politik["Konservatif", "Bush"]) /
(tabel_politik["Liberal", "Bush"] * tabel_politik["Konservatif", "Dukakis"])
OR
## [1] 10.78597
Berdasarkan hasil analisis log-linear:
Model independen (tanpa interaksi) memiliki nilai Deviance sebesar 170.4841 dengan p-value = 0, yang berarti model ini tidak cocok secara signifikan terhadap data.
Sebaliknya, model saturated (yang mencakup interaksi antara pandangan politik dan pilihan presiden) menunjukkan deviance = 0, artinya cocok sempurna dengan data.
Karena nilai p-value = 0 untuk selisih deviance, maka interaksi antara pandangan politik dan pilihan presiden signifikan. Dengan kata lain, pilihan kandidat presiden sangat bergantung pada orientasi politik responden.
Terdapat asosiasi yang sangat kuat antara orientasi politik seseorang dan pilihan mereka dalam pemilu presiden 1988. Interaksi ini tidak hanya signifikan secara statistik, tetapi juga besar secara praktis, seperti yang terlihat dari nilai odds ratio yang tinggi. Ini menunjukkan pentingnya mempertimbangkan efek interaksi dalam analisis data kategorik politik.
Agresti, A. (2002). Categorical Data Analysis. John Wiley & Sons.
Kleinbaum, D. G., & Klein, M. (2010). Logistic Regression: A Self-Learning Text. Springer.
Hosmer, D. W., Lemeshow, S., & Sturdivant, R. X. (2013). Applied Logistic Regression. John Wiley & Sons.
Sudjana, N. (2005). Metoda Statistika. Tarsito.
Ghozali, I. (2016). Aplikasi Analisis Multivariat dengan Program IBM SPSS 23. Badan Penerbit Universitas Diponegoro.
Utts, J. M., & Heckard, R. F. (2015). Mind on Statistics. Cengage Learning.
International Diabetes Federation (IDF). (2021). *IDF Diabetes Atlas (10th ed.)
Kahn, B. B., & Flier, J. S. (2006). Obesity and insulin resistance. Journal of Clinical Investigation, 106(4), 473–481.
Kementerian Kesehatan Republik Indonesia. (2018). Laporan Nasional Riskesdas 2018.
Menke, A., Casagrande, S., Geiss, L., & Cowie, C. C. (2015). Prevalence of and trends in diabetes among adults in the United States, 1988-2012. JAMA, 314(10), 1021–1029.
Anderson, C. (n.d.). Log-linear Models for Contingency
Tables. University of Illinois.
Diakses dari: https://education.illinois.edu/docs/default-source/carolyn-anderson/edpsy589/lectures/6-loglinear/6_loglinear_models_beamer_post.pdf