- Bab I Pendahuluan
- Bab II Metode dalam Analisis Data Kategori
- Bab III Distribusi Probabilitas dalam Data Kategori
- Bab IV Desain Sampling dalam Analisis Data Kategori
- Bab V Tabel Kontingensi 2 × 2
- Bab VI Inferensi Tabel Kontingensi Dua Arah
- Bab 7 Tabel Kontingensi Tiga Arah
- Bab 8 Generalized Linear Model (GLM)
- Bab IX Inferensi GLM
- Bab X Regresi Logistik dengan Prediktor Nominal, Ordinal, dan Rasio
- Bab XI
Pemilihan Model Regresi Logistik dan Evaluasi
- 11.1 Membangun Model Regresi Logistik: Pendekatan Confirmatory dan Exploratory
- 11.2 Pemilihan Model
- 11.3 Evaluasi Model: ROC dan AUC
- 11.4 Pseudo R-Squared
- 11.5 Tabel Klasifikasi dan Evaluasi
- 11.6 Metode Perbandingan Model dalam Regresi Logistik
- 11.7 Likelihood-Ratio Test
- 11.8 Prinsip Parsimony
- 11.9 Evaluasi Tabel Klasifikasi dan Akurasi Model
- 11.10 ROC (Receiver Operating Characteristic)
- 11.11 Precision-Recall Curve (PR Curve)
- 11.12 Pseudo R-squared pada Regresi Logistik
- Bab XII Distribusi Multinomial
- Bab XIII Regresi Logistik Ordinal
- Bab XIV Model Log Linear Dua Arah
- Bab XV Model Log Linear Tiga Arah
- Bab XVI Uji Model Interaksi Tiga Arah (Saturated vs Homogenous)
- Bab XVII Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on X)
- Bab XVIII Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on Y)
- Bab XIX Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on Z)
- Bab XX Pemilihan Model Terbaik
- Studi Kasus 1
- Studi Kasus 2
- Referensi
Bab I Pendahuluan
1.1 Pengertian
Analisis data kategori adalah metode statistik yang digunakan untuk menganalisis data yang diklasifikasikan ke dalam kategori-kategori, seperti jenis kelamin, warna, status, atau kelas. Berbeda dari data numerik, data ini tidak memiliki nilai angka yang memiliki makna matematis, melainkan menunjukkan kelompok atau label (Agresti, 2018).
1.2 Tujuan
Tujuan utama dari analisis data kategori adalah untuk mengeksplorasi dan menguji hubungan antara variabel-variabel kategorik, membuat prediksi, serta memahami pola distribusi dalam populasi berdasarkan data klasifikasi. Ini juga membantu pengambilan keputusan berbasis data pada berbagai bidang (Agresti, 2002).
1.3 Definisi dan Ruang Lingkup
Analisis ini mencakup metode-metode seperti tabel kontingensi, uji chi-square untuk independensi, correspondence analysis untuk visualisasi, hingga pendekatan machine learning seperti decision tree dan random forest. Fokusnya adalah pada variabel-variabel yang dikodekan sebagai kategori (McHugh, 2013).
1.4 Manfaat dalam Berbagai Bidang
Analisis data kategori memiliki penerapan luas di berbagai bidang karena banyak variabel dalam kehidupan nyata bersifat kategorik. Berikut ini beberapa bidang yang mendapat manfaat besar:
1.4.1 Kesehatan dan Epidemiologi
Berguna untuk mengidentifikasi hubungan antara faktor risiko (seperti merokok, jenis kelamin, usia) dengan kejadian penyakit tertentu. Contohnya, uji chi-square sering digunakan untuk menganalisis prevalensi penyakit berdasarkan jenis kelamin atau kelompok usia (McHugh, 2013).
1.4.2 Ilmu Sosial dan Politik
Digunakan untuk menganalisis hubungan antara atribut sosial (seperti umur, jenis kelamin, status sosial) dan sikap atau perilaku, misalnya preferensi politik atau tingkat kepercayaan terhadap pemerintah (Agresti, 2002).
1.4.3 Pemasaran dan Bisnis
Perusahaan menganalisis preferensi pelanggan berdasarkan kategori demografis untuk menyusun strategi penjualan dan promosi yang tepat sasaran. Misalnya, klasifikasi pelanggan berdasarkan tingkat pendapatan atau lokasi untuk personalisasi iklan (Agresti, 2018).
1.4.4 Teknologi dan Data Science:
Model seperti decision tree dan random forest digunakan untuk klasifikasi dalam sistem rekomendasi, pendeteksian spam, atau pengenalan pola perilaku pengguna berdasarkan atribut kategorikal seperti jenis perangkat atau waktu akses (Breiman, 2001).
1.4.5 Pendidikan
Digunakan untuk menganalisis hasil belajar berdasarkan kategori seperti latar belakang pendidikan, jenis kelamin, dan metode pengajaran. Hal ini membantu dalam merancang kebijakan pendidikan yang lebih inklusif (Agresti, 2018).
1.4.6 Psikologi
Kuesioner psikologis yang menghasilkan data kategorik (seperti ya/tidak, skala likert) dapat dianalisis untuk mengetahui hubungan antara variabel seperti tingkat stres dan faktor lingkungan (Greenacre, 2007)
1.4.7 Hukum dan Kriminologi
Digunakan untuk menganalisis data kriminal berdasarkan jenis kejahatan, lokasi geografis, atau profil pelaku. Ini mendukung kebijakan keamanan berbasis data (Agresti, 2002)
Bab II Metode dalam Analisis Data Kategori
2.1 Tabel Kontingensi
Tabel kontingensi menyajikan frekuensi atau jumlah observasi untuk kombinasi dua atau lebih variabel kategorik. Ini adalah dasar analisis hubungan antar kategori dan memberikan gambaran awal mengenai asosiasi atau ketidakteraturan distribusi data (Agresti, 2002).
2.2 Uji Chi-Square
Uji chi-square digunakan untuk menguji apakah ada hubungan yang signifikan antara dua variabel kategorik. Hipotesis nol biasanya menyatakan bahwa kedua variabel tidak berhubungan (independen), sedangkan hipotesis alternatif menyatakan adanya hubungan (McHugh, 2013).
2.3 Analisis Correspondence (CA)
Correspondence Analysis (CA) adalah metode eksploratori untuk menganalisis dan memvisualisasikan hubungan antara dua variabel kategorik dalam bentuk peta dua dimensi. CA sering digunakan untuk menyederhanakan interpretasi data yang kompleks dalam tabel kontingensi besar (Greenacre, 2007).
2.4 Decision Tree dan Random Forest
2.4.1 Decision Tree
Algoritma pembelajaran terarah (supervised learning) yang membentuk model prediksi berbasis aturan keputusan dari data kategorik atau numerik.
2.4.2 Random Forest
Ensembel dari banyak decision tree yang digabungkan untuk meningkatkan akurasi dan mengurangi overfitting. Cocok untuk klasifikasi dan regresi pada data kompleks (Breiman, 2001).
Bab III Distribusi Probabilitas dalam Data Kategori
Distribusi probabilitas dalam data kategori digunakan untuk memodelkan kemungkinan hasil dari percobaan yang menghasilkan data diskret. Distribusi-distribusi ini adalah fondasi teori dalam statistik kategorikal.
3.1 Distribusi Bernoulli
Distribusi Bernoulli digunakan untuk percobaan dengan dua kemungkinan
hasil: “sukses” (1) atau “gagal” (0). Misalnya, apakah seseorang
terdiagnosis hipertensi (ya/tidak). Distribusi ini memiliki satu
parameter: probabilitas sukses𝑝(Agresti, 2018). \[
P(X = x) = p^x (1 - p)^{1 - x}, \quad x \in \{0, 1\}
\] - \(p\): probabilitas
sukses
- \(1 - p\): probabilitas gagal
Contoh: Seorang pasien minum obat (1) atau tidak (0).
3.2 Distribusi Binomial
Distribusi Binomial adalah perluasan dari distribusi Bernoulli ke
𝑛pengulangan (Hosmer, Lemeshow, & Sturdivant, 2013).
\[
P(X = x) = \binom{n}{x} p^x (1 - p)^{n - x}, \quad x = 0, 1, ..., n
\] - \(n\): jumlah
percobaan
- \(x\): jumlah sukses
- \(p\): probabilitas sukses
Contoh: Dari 10 pasien, berapa yang rutin minum obat jika probabilitas kepatuhan adalah 0,7?
3.3 Distribusi Multinomial
Distribusi Multinomial adalah perluasan dari Binomial ketika hasil tidak hanya dua kategori, melainkan lebih dari dua (Agresti, 2002). \[ P(X_1 = x_1, \dots, X_k = x_k) = \frac{n!}{x_1! \cdots x_k!} p_1^{x_1} \cdots p_k^{x_k} \]
Dengan syarat: \[
\sum_{i=1}^{k} x_i = n, \quad \sum_{i=1}^{k} p_i = 1
\] - \(x_i\): jumlah kejadian
pada kategori ke-\(i\)
- \(p_i\): probabilitas kategori
ke-\(i\)
- \(n\): total pengamatan
Contoh: Respon pasien terhadap terapi: membaik, tetap, atau memburuk.
3.4 Distribusi Poisson
Distribusi Poisson digunakan untuk menghitung jumlah kejadian dalam
interval waktu atau ruang tertentu (Dobson & Barnett, 2008).
\[
P(X = x) = \frac{e^{-\lambda} \lambda^x}{x!}, \quad x = 0, 1, 2, ...
\] - \(\lambda\): rata-rata
kejadian dalam interval
- \(e\): konstanta Euler (\(\approx 2.718\))
Contoh: Jumlah pasien hipertensi yang datang ke puskesmas per hari.
Bab IV Desain Sampling dalam Analisis Data Kategori
Desain sampling dalam analisis data kategori sangat penting untuk validitas hasil dan jenis inferensi yang bisa ditarik dari data. Terdapat dua pendekatan utama yaitu prospektif dan retrospektif.
4.1 Prospective Sampling
Desain prospektif dilakukan dengan mengikuti subjek dari waktu ke waktu ke depan, biasanya dimulai dari kondisi awal (misal paparan) dan mencatat hasilnya di masa depan.
4.1.1 Eksperimen
Eksperimen adalah desain sampling di mana subjek secara acak dibagi
ke dalam kelompok perlakuan dan kontrol.
Tujuan utamanya adalah untuk menguji hubungan kausal antara perlakuan
dan outcome.
Contoh: Uji coba acak terkontrol (Randomized Controlled Trial/RCT) untuk menguji efektivitas terapi hipertensi baru.
4.1.2 Studi Kohort
Dalam studi kohort, subjek diklasifikasikan berdasarkan status paparan terhadap suatu faktor risiko, dan kemudian diikuti ke depan untuk melihat apakah mereka mengalami outcome.
Contoh: Studi terhadap perokok dan non-perokok untuk mengamati kejadian hipertensi selama 10 tahun.
4.2 Retrospective Sampling
Sampling retrospektif dilakukan dengan mengumpulkan data dari kejadian yang sudah terjadi (misalnya, dari rekam medis).
4.2.1 Studi Kasus-Kontrol
Dalam studi kasus-kontrol, subjek dipilih berdasarkan outcome (misal, pasien yang sudah mengalami hipertensi dan yang tidak), lalu dicari informasi ke belakang mengenai paparan mereka.
Contoh: Studi yang melihat riwayat konsumsi garam pada pasien dengan dan tanpa hipertensi.
4.2.2 Studi Kohort Retrospektif
Desain ini mirip dengan kohort, tetapi semua data dikumpulkan dari dokumen atau rekam masa lalu.
Contoh: Meneliti hubungan antara aktivitas fisik dan kejadian stroke berdasarkan data rumah sakit 5 tahun terakhir.
4.3 Tabel Perbandingan Desain Sampling
Berikut adalah ringkasan perbedaan antar desain:
## Warning: package 'knitr' was built under R version 4.3.3| Jenis.Desain | Waktu.Pengamatan | Kelebihan | Kekurangan |
|---|---|---|---|
| Eksperimen (RCT) | Prospektif | Kontrol tinggi terhadap variabel, cocok untuk inferensi kausal | Biaya tinggi, etika dan waktu pelaksanaan lama |
| Studi Kohort Prospektif | Prospektif | Dapat mengukur insidensi dan risiko relatif secara langsung | Memerlukan waktu lama dan rawan loss to follow-up |
| Studi Kasus-Kontrol | Retrospektif | Efisien untuk penyakit langka dan murah | Tidak bisa mengukur insidensi, rentan bias recall |
| Studi Kohort Retrospektif | Retrospektif | Cepat dan murah karena menggunakan data historis | Rentan terhadap bias dan ketidaklengkapan data |
Bab V Tabel Kontingensi 2 × 2
Tabel kontingensi 2 × 2 digunakan untuk menampilkan frekuensi dari dua variabel kategorikal yang masing-masing memiliki dua kategori. Tabel ini sangat umum digunakan dalam studi epidemiologi, uji klinis, dan analisis hubungan antara dua variabel diskrit. Berikut tabel kontingensi 2x2:
| Jenis Kelamin | Yes | No | Total |
|---|---|---|---|
| Female | 371 | 74 | 445 |
| Male | 250 | 71 | 321 |
| Total | 621 | 145 | 766 |
5.1 Distribusi Peluang dalam Tabel Kontingensi 2 × 2
5.1.1 Peluang Bersama (Joint Probability)
Peluang bersama adalah probabilitas bahwa dua kejadian terjadi secara bersamaan (Agresti, 2018).
\[ P(A \cap B) = \frac{n_{AB}}{n} \] - \(n_{AB}\): jumlah observasi dengan kejadian A dan B - \(n\): total observasi
5.1.2 Peluang Marginal (Marginal Probability)
Peluang marginal adalah probabilitas terjadinya satu kejadian, tanpa memperhatikan kejadian lainnya (Agresti, 2018).
\[ P(A) = \frac{\sum_B n_{AB}}{n}, \quad P(B) = \frac{\sum_A n_{AB}}{n} \]
5.1.3 Peluang Bersyarat (Conditional Probability)
Peluang bersyarat adalah probabilitas kejadian A, dengan syarat kejadian B sudah terjadi (McClave & Sincich, 2014).
\[ P(A | B) = \frac{P(A \cap B)}{P(B)} = \frac{n_{AB}}{n_{.B}} \] - \(n_{.B}\): jumlah kasus dengan B - (A | B): peluang A jika diketahui B
Perhitungan manual
Peluang Bersama (Joint Probability)
Peluang bersama dihitung dengan membagi frekuensi sel dengan total observasi.
\[ \begin{aligned} P(\text{Yes, Female}) &= \frac{371}{766} \approx 0.4842 \\ P(\text{No, Female}) &= \frac{74}{766} \approx 0.0966 \\ P(\text{Yes, Male}) &= \frac{250}{766} \approx 0.3265 \\ P(\text{No, Male}) &= \frac{71}{766} \approx 0.0927 \\ \end{aligned} \]
Peluang Marginal
\[ \begin{aligned} P(\text{Female}) &= \frac{445}{766} \approx 0.5811 \\ P(\text{Male}) &= \frac{321}{766} \approx 0.4189 \\ P(\text{Yes}) &= \frac{621}{766} \approx 0.8114 \\ P(\text{No}) &= \frac{145}{766} \approx 0.1892 \\ \end{aligned} \]
Peluang Bersyarat
\[ \begin{aligned} P(\text{Yes | Female}) &= \frac{371}{445} \approx 0.8337 \\ P(\text{No | Female}) &= \frac{74}{445} \approx 0.1663 \\ P(\text{Yes | Male}) &= \frac{250}{321} \approx 0.7782 \\ P(\text{No | Male}) &= \frac{71}{321} \approx 0.2218 \\ \end{aligned} \]
Perhitungan dengan R
# Buat data observasi
data <- matrix(c(371, 74, 250, 71), nrow = 2, byrow = TRUE)
colnames(data) <- c("Yes", "No")
rownames(data) <- c("Female", "Male")
n <- sum(data)
# Peluang Bersama
P_joint <- data / n
P_joint## Yes No
## Female 0.4843342 0.09660574
## Male 0.3263708 0.09268930# Peluang Marginal
P_marginal_rows <- rowSums(data) / n
P_marginal_cols <- colSums(data) / n
P_marginal_rows## Female Male
## 0.5809399 0.4190601P_marginal_cols## Yes No
## 0.810705 0.189295# Peluang Bersyarat
P_conditional <- prop.table(data, margin = 1)
P_conditional## Yes No
## Female 0.8337079 0.1662921
## Male 0.7788162 0.2211838Interpretasi:
Peluang Bersama: peluang seseorang yang dipilih secara acak adalah perempuan dan menjawab “Yes” adalah sekitar 48.4%. Peluang seseorang adalah laki-laki dan menjawab “No” adalah sekitar 9.3%.
Peluang Marginal: mayoritas responden adalah perempuan. Sebagian besar responden menjawab “Yes”
Peluang Bersyarat: Dari seluruh perempuan, sekitar 83.4% menjawab “Yes”. Dari seluruh laki-laki, sekitar 77.9% menjawab “Yes”. Proporsi Yes lebih tinggi pada perempuan dibandingkan laki-laki.
5.2 Ukuran Asosiasi dalam Data Kategori 2 × 2
Dalam analisis statistik, tabel kontingensi 2×2 digunakan untuk mengevaluasi hubungan antara dua variabel kategorik. Tabel ini menyajikan frekuensi gabungan dari dua variabel yang masing-masing memiliki dua kategori, sehingga memungkinkan kita untuk mengidentifikasi ada tidaknya keterkaitan statistik di antara keduanya.Pendekatan ini umum digunakan di berbagai disiplin ilmu, seperti:
Epidemiologi: Untuk mengevaluasi kaitan antara perilaku merokok dan kejadian kanker paru-paru.
Uji Klinis: Dalam menilai efektivitas suatu pengobatan terhadap kondisi medis tertentu.
Ilmu Sosial: Untuk meneliti hubungan antara jenjang pendidikan dengan tingkat pekerjaan.
5.2.1 Risk Difference (RD)
RRisk Diference (RD) atau Perbedaan Risiko menunjukkan selisih risiko antara kelompok terpapar dan tidak terpapar (Rothman et al., 2008).
\[ RD = \frac{a}{a + b} - \frac{c}{c + d} \]
Jika RD > 0, maka risiko kejadian lebih tinggi di Grup 1 dibandingkan Grup 2.
Jika RD < 0, maka risiko kejadian lebih rendah di Grup 1 dibanding
Jika RD = 0, maka tidak ada perbedaan risiko antara dua kelompok.
5.2.2 Relative Risk (RR)
Relative risk menunjukkan perbandingan risiko antara kelompok terpapar dan tidak terpapar (Kleinbaum et al., 2007).
\[ RR = \frac{ \frac{a}{a + b} }{ \frac{c}{c + d} } \]
Jika RR > 1, maka kejadian lebih sering terjadi di Grup 1 dibandingkan Grup 2.
Jika RR < 1, maka kejadian lebih jarang terjadi di Grup 1 dibanding
Jika RR = 1, maka tidak ada perbedaan risiko antara dua kelompok.
5.2.3 Odds Ratio (OR)
Odds ratio menunjukkan perbandingan peluang terjadinya kejadian pada dua kelompok (Agresti, 2002).
\[ OR = \frac{a \cdot d}{b \cdot c} \]
Jika OR > 1, maka peluang kejadian lebih besar di Grup 1 dibandingkan Grup 2.
Jika OR < 1, maka peluang kejadian lebih kecil di Grup 1 dibandingkan Grup 2.
Jika OR 1, maka tidak ada perbedaan peluang kejadian antara dua kelompok
| Ukuran | Definisi | Desain.Sampling.yang.Cocok | Interpretasi |
|---|---|---|---|
| Risk Difference (RD) | Selisih probabilitas kejadian antara dua kelompok | Studi kohort atau eksperimen acak | Menunjukkan tambahan atau pengurangan risiko absolut |
| Relative Risk (RR) | Perbandingan risiko kejadian di dua kelompok | Studi kohort atau eksperimen klinis | Mengukur berapa kali lebih besar risiko di satu kelompok dibandingkan kelompok lainnya |
| Odds Ratio (OR) | Perbandingan odds antara dua kelompok | Studi kasus-kontrol atau studi observasional | Mengukur hubungan antara variabel eksposur dan kejadian dalam desain studi kasus-kontrol |
Perhitungan manual
Hitung Risiko
Risiko pada kelompok Female (terpapar):
\[ \text{Risk}_{Female} = \frac{a}{a + b} = \frac{371}{371 + 74} = \frac{371}{445} \approx 0.8337 \]
Risiko pada kelompok Male (tidak terpapar):
\[ \text{Risk}_{Male} = \frac{c}{c + d} = \frac{250}{250 + 71} = \frac{250}{321} \approx 0.7788 \]
Risk Difference (RD)
\[ RD = \text{Risk}_{Female} - \text{Risk}_{Male} = 0.8337 - 0.7788 = 0.0549 \] Interpretasi:
Risiko percaya pada kehidupan setelah mati lebih tinggi sekitar 5.49% pada perempuan dibanding laki-laki.
Relative Risk (RR)
\[ RR = \frac{\text{Risk}_{Female}}{\text{Risk}_{Male}} = \frac{0.8337}{0.7788} \approx 1.0705 \] Interpretasi:
Perempuan memiliki peluang sekitar 1.07 kali lebih besar untuk percaya pada kehidupan setelah mati dibanding laki-laki.
Odds Ratio (OR)
Odds pada Female:
\[ \text{Odds}_{Female} = \frac{a}{b} = \frac{371}{74} \approx 5.0135 \]
Odds pada Male:
\[ \text{Odds}_{Male} = \frac{c}{d} = \frac{250}{71} \approx 3.5211 \]
Maka:
\[ OR = \frac{a \cdot d}{b \cdot c} = \frac{371 \cdot 71}{74 \cdot 250} \approx \frac{26,341}{18,500} \approx 1.4238 \] Interpretasi:
Perempuan memiliki odds sekitar 1.42 kali lebih tinggi untuk percaya pada kehidupan setelah mati dibanding laki-laki.
Perhitungan dengan R
## [1] 0.05489167## [1] 1.070481## [1] 1.423838Interpretasi:
Perempuan memiliki odds sekitar 1.42 kali lebih tinggi untuk percaya pada kehidupan setelah mati dibanding laki-laki.
Bab VI Inferensi Tabel Kontingensi Dua Arah
Inferensi dalam statistik mengacu pada proses pengambilan kesimpulan mengenai populasi berdasarkan sampel data. Dalam konteks tabel kontingensi dua arah, inferensi digunakan untuk menganalisis hubungan antara dua variabel kategorikal yang disusun dalam tabel kontingensi.
Tabel kontingensi adalah tabel yang menyajikan distribusi frekuensi dari dua variabel kategorikal dalam bentuk matriks. Tujuan utama dari inferensi pada tabel ini adalah untuk memahami hubungan antara variabel-variabel tersebut melalui estimasi dan pengujian hipotesis. Inferensi dalam tabel kontingensi dua arah dapat dibagi menjadi dua kategori utama, yaitu Estimasi dan Pengujian
6.1 Estimasi
Estimasi bertujuan untuk memperkirakan parameter populasi berdasarkan data sampel. Estimasi dibagi menjadi:
6.1.1 Estimasi Titik
Estimasi titik adalah nilai tunggal yang diperkirakan sebagai parameter populasi. Dalam tabel kontingensi, proporsi adalah contoh umum (Agresti, 2018). \[ \hat{p} = \frac{x}{n} \]
- \(x\): jumlah kejadian sukses
- \(n\): jumlah total observasi
6.1.2 Estimasi Interval
Estimasi interval memberikan rentang nilai untuk parameter populasi, biasanya dalam bentuk interval kepercayaan/confidence interval (Agresti, 2018). Interval kepercayaan 95% untuk proporsi:
\[ CI = \hat{p} \pm z_{\alpha/2} \cdot \sqrt{ \frac{\hat{p}(1 - \hat{p})}{n} } \]
- \(z_{\alpha/2}\): nilai kritis distribusi normal (biasanya 1.96 untuk 95% CI)
- \(hat{p}\): estimasi titik proporsi
- \(n\): ukuran sampel
Perhitungan dengan R
## [1] 0.05489167## [1] -0.002189698 0.111973029## [1] 1.070481## [1] 0.9965561 1.1498894## [1] 1.423838## [1] 0.9904462 2.0468696Interpretasi:
Meskipun nilai RD, RR, dan OR menunjukkan arah efek yang mengarah pada risiko/peluang lebih tinggi pada perempuan, tidak satupun dari ketiganya signifikan secara statistik pada level 95%, karena RD CI mencakup 0 dan RR & OR CI mencakup 1
Dengan kata lain, tidak ada cukup bukti bahwa jenis kelamin berhubungan signifikan dengan kejadian yang diamati (misalnya, kejadian penyakit, partisipasi, atau outcome lain tergantung konteks datanya).
6.2 Uji Hipotesis
6.2.1 Uji Proporsi
Formulasi Uji Proporsi
Untuk menguji hipotesis bahwa tidak ada perbedaan proporsi antara dua kelompok, kita menggunakan uji z dua proporsi, dengan hipotesis:
- Hipotesis Nol (\(H_0\)): Tidak ada perbedaan proporsi antara dua kelompok, yaitu \(p_1 = p_2\)
- Hipotesis Alternatif (\(H_1\)): Terdapat perbedaan proporsi antara dua kelompok, yaitu \(p_1 \ne p_2\)
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 (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.
Perhitungan Manual
Proporsi Sampel
- Proporsi untuk kelompok Female:
\[ \hat{p}_1 = \frac{371}{445} = 0.8337 \]
- Proporsi untuk kelompok Male:
\[ \hat{p}_2 = \frac{250}{321} = 0.7788 \]
Proporsi Gabungan
Gabungkan jumlah “Yes” dari kedua kelompok:
\[ \hat{p} = \frac{371 + 250}{445 + 321} = \frac{621}{766} = 0.8111 \]
Statistik Uji Z
\[ Z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{\hat{p}(1 - \hat{p}) \left( \frac{1}{n_1} + \frac{1}{n_2} \right)}} \]
\[ Z = \frac{0.8337 - 0.7788}{\sqrt{0.8111(1 - 0.8111)\left(\frac{1}{445} + \frac{1}{321}\right)}} \]
\[ Z = \frac{0.0549}{\sqrt{0.8111 \cdot 0.1889 \left(0.002247 + 0.003115\right)}} \]
\[ Z = \frac{0.0549}{\sqrt{0.8111 \cdot 0.1889 \cdot 0.005362}} = \frac{0.0549}{\sqrt{0.000820}} = \frac{0.0549}{0.0286} \approx 1.92 \]
Perhitungan dengan R
## [1] 0.8337079## [1] 0.7788162## [1] 0.810705## [1] 1.913478Interpretasi:
Meskipun proporsi respon “Yes” pada perempuan (83.37%) sedikit lebih tinggi dibandingkan laki-laki (77.88%), perbedaan ini tidak cukup signifikan secara statistik pada tingkat signifikansi 5%.
6.2.2 Uji Asosiasi
Untuk setiap uji asosiasi, hipotesis yang diuji adalah:
- Hipotesis Nol (\(H_0\)): Tidak ada asosiasi antara dua variabel.
- Hipotesis Alternatif (\(H_1\)): Terdapat asosiasi antara dua variabel.
6.2.2.1 Risk Difference (RD)
Risk Diference mengukur perbedaan absolut dalam probabilitas kejadian antara dua kelompok:
\[ RD = \frac{n_{11}}{n_1} - \frac{n_{21}}{n_2} \]
Standard Error:
\[ SE(RD) = \sqrt{ \frac{\hat{p}_1(1 - \hat{p}_1)}{n_1} + \frac{\hat{p}_2(1 - \hat{p}_2)}{n_2} } \]
Z-statistik:
\[ Z_{RD} = \frac{RD}{SE(RD)} \]
6.2.2.2 Relative Risk (RR)
Relative Risk membandingkan kemungkinan kejadian antara dua kelompok:
\[ RR = \frac{n_{11} / n_1}{n_{21} / n_2} \]
SE log(RR):
\[ SE(\ln RR) = \sqrt{ \frac{1}{n_{11}} - \frac{1}{n_1} + \frac{1}{n_{21}} - \frac{1}{n_2} } \]
Z-statistik:
\[ Z_{RR} = \frac{\ln RR}{SE(\ln RR)} \]
6.2.2.3 Odds Ratio (OR)
Odds Ratio membandingkan peluang kejadian antara dua kelompok:
\[ OR = \frac{n_{11} \cdot n_{22}}{n_{12} \cdot n_{21}} \]
SE log(OR):
\[ SE(\ln OR) = \sqrt{ \frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}} } \]
Z-statistik:
\[ Z_{OR} = \frac{\ln OR}{SE(\ln OR)} \]
Perhitungan dengan R
## $Risk_Difference
## [1] 0.05489167
##
## $SE_RD
## [1] 0.02912368
##
## $Z_RD
## [1] 1.884778## $Relative_Risk
## [1] 1.070481
##
## $SE_log_RR
## [1] 0.03650975
##
## $Z_RR
## [1] 1.865474## $Odds_Ratio
## [1] 1.423838
##
## $SE_log_OR
## [1] 0.1851849
##
## $Z_OR
## [1] 1.908125Interpretasi:
Berdasarkan hasil analisis, diperoleh bahwa Risk Difference (RD) antara kelompok perempuan dan laki-laki adalah sebesar 0.0549 dengan standard error (SE) sebesar 0.0291 dan nilai statistik Z = 1.88, yang menunjukkan perbedaan risiko tidak signifikan secara statistik pada tingkat signifikansi 5% (karena Z < 1.96). Relative Risk (RR) sebesar 1.0705 menunjukkan bahwa perempuan memiliki risiko 7% lebih tinggi dibandingkan laki-laki, namun dengan SE log(RR) = 0.0365 dan Z = 1.87, perbedaan ini juga tidak signifikan. Sementara itu, Odds Ratio (OR) sebesar 1.4238 menunjukkan bahwa odds perempuan untuk mengalami kejadian sekitar 42% lebih tinggi dibandingkan laki-laki, tetapi dengan SE log(OR) = 0.1852 dan Z = 1.91, hasil ini juga belum mencapai batas signifikansi statistik (Z < 1.96). Dengan demikian, meskipun ada kecenderungan perbedaan antara kelompok, tidak terdapat bukti statistik yang cukup kuat untuk menyimpulkan adanya perbedaan yang signifikan.
6.2.3 Uji Independensi
Uji independensi digunakan untuk menentukan apakah ada hubungan statistik antara dua variabel kategorikal.
6.2.3.1 Uji Chi-Square
Uji Chi-Square digunakan untuk menguji apakah ada hubungan antara dua variabel kategorikal.
Rumus Chi-Square: \[ \chi^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:
- \(R_i\) = total baris ke-i
- \(C_j\) = total kolom ke-j
- \(N\) = total sampel
Perhitungan dengan R
##
## Pearson's Chi-squared test
##
## data: data
## X-squared = 3.6614, df = 1, p-value = 0.05569Interpretasi:
Hasil uji Chi-Square menunjukkan nilai statistik χ² = 3.6614 dengan derajat bebas (df) sebesar 1 dan p-value = 0.05569.
Karena p-value (0.05569) > 0.05, maka kita tidak menolak H₀ (hipotesis nol). Artinya, tidak terdapat bukti yang cukup untuk menyatakan adanya hubungan yang signifikan secara statistik antara jenis kelamin dan kejadian (ya/tidak) pada tingkat signifikansi 5%.
Namun, perlu dicatat bahwa p-value ini mendekati 0.05, sehingga bisa dikatakan ada indikasi lemah adanya hubungan, tetapi belum cukup kuat secara statistik untuk dinyatakan signifikan.
6.2.3.2 Partisi Chi-Square
Partisi Chi-Square digunakan untuk mengidentifkasi kategori mana dalam tabel kontingensi yang bertanggung jawab atas hubungan yang signifkan. Jika uji Chi-Square pada tabel kontingensi I×J signifkan, maka partisi Chi-Square memungkinkan kita untuk menguraikan efek hubungan dalam subkelompok yang lebih kecil.
Menurut Lancaster dan Irwin, statistik uji chi-square dapat dipecah menjadi komponen-komponen sebanyak derajat bebas dari tabel kontingensi, memungkinkan identifkasi kategori yang berkontribusi pada hubungan yang signifkan.
Hal ini juga dikenal dengan Simpson’s Paradox: Tren yang muncul dalam beberapa kelompok data dapat menghilang atau bahkan berbalik arah ketika data tersebut digabungkan.
| Gender | Yes | Undecided | No |
|---|---|---|---|
| White | 435 | 58 | 89 |
| Black | 275 | 50 | 84 |
| Total | 710 | 108 | 173 |
Hipotesis yang diuji adalah:
Hipotesis Nol: Tidak ada hubungan antara variabel Gender dan Identifkasi Partai Politik.
Hipotesis Alternatif: Ada hubungan antara variabel Gender dan Identifkasi Partai Politik.
Perhitungan dengan R
# Data Observasi
data_matrix <- matrix(c(435, 58, 89, 275, 50, 84), nrow = 2, byrow = TRUE)
colnames(data_matrix) <- c("Yes", "Undecided", "No")
rownames(data_matrix) <- c("White", "Black")
# Uji Chi-Square
chi_test <- chisq.test(data_matrix)
# Hasil
list(
Chi_Square = chi_test$statistic,
P_Value = chi_test$p.value,
Decision = ifelse(chi_test$p.value < 0.05, "Tolak H0", "Gagal Tolak H0")
)## $Chi_Square
## X-squared
## 6.799858
##
## $P_Value
## [1] 0.03337564
##
## $Decision
## [1] "Tolak H0"# Data Partisi 1: Yes vs No
data_partisi1 <- matrix(c(435, 89, 275, 84), nrow = 2, byrow = TRUE)
colnames(data_partisi1) <- c("Yes", "No")
rownames(data_partisi1) <- c("White", "Black")
# Uji Chi-Square Partisi 1
chi_test1 <- chisq.test(data_partisi1)
# Data Partisi 2: (Yes + No) vs Undecided
yes_no_white <- 435 + 89
yes_no_black <- 275 + 84
data_partisi2 <- matrix(c(yes_no_white, 58, yes_no_black, 50), nrow = 2, byrow = TRUE)
colnames(data_partisi2) <- c("Yes+No", "Undecided")
rownames(data_partisi2) <- c("White", "Black")
# Uji Chi-Square Partisi 2
chi_test2 <- chisq.test(data_partisi2)
# Hasil
list(
Chi_Square_Partisi1 = list(
Statistic = chi_test1$statistic,
P_Value = chi_test1$p.value,
Decision = ifelse(chi_test1$p.value < 0.05, "Tolak H0", "Gagal Tolak H0")
),
Chi_Square_Partisi2 = list(
Statistic = chi_test2$statistic,
P_Value = chi_test2$p.value,
Decision = ifelse(chi_test2$p.value < 0.05, "Tolak H0", "Gagal Tolak H0")
)
)## $Chi_Square_Partisi1
## $Chi_Square_Partisi1$Statistic
## X-squared
## 5.163027
##
## $Chi_Square_Partisi1$P_Value
## [1] 0.02307266
##
## $Chi_Square_Partisi1$Decision
## [1] "Tolak H0"
##
##
## $Chi_Square_Partisi2
## $Chi_Square_Partisi2$Statistic
## X-squared
## 1.040704
##
## $Chi_Square_Partisi2$P_Value
## [1] 0.3076577
##
## $Chi_Square_Partisi2$Decision
## [1] "Gagal Tolak H0"Interpretasi:
Uji Chi-Square ini digunakan untuk mengetahui apakah ada hubungan antara dua variabel kategorik, yaitu Gender (White/Black) dan Respons (Yes/Undecided/No). Dengan P-Value = 0.033, yang lebih kecil dari 0.05, kita menolak H0 pada tingkat signifikansi 5%.
Terdapat hubungan yang signifikan secara statistik antara ras (White vs. Black) dan pilihan jawaban (Yes, Undecided, No). Artinya, distribusi jawaban tidak sama antara responden White dan Black — misalnya, satu kelompok mungkin lebih cenderung menjawab Yes dibanding yang lain.
Terdapat perbedaan yang signifikan secara statistik dalam pola jawaban Yes dan No antara responden berkulit putih dan hitam. Ini berarti bahwa distribusi pilihan antara “Yes” dan “No” tidak merata antar ras, atau dalam kata lain, ras berkaitan dengan jawaban yang diberikan (antara setuju dan tidak setuju).
Tidak ada bukti yang cukup untuk menyatakan bahwa terdapat perbedaan signifikan antara kelompok ras dalam hal memilih “Undecided” dibandingkan “Yes/No”. Ini menunjukkan bahwa keraguan (Undecided) terhadap isu tersebut tidak terlalu dipengaruhi oleh ras.
6.2.3.2 Uji Likelihood Ratio (G²)
Uji Likelihood Ratio (G²) adalah alternatif dari uji chi-square yang digunakan untuk menguji hipotesis independensi dalam tabel kontingensi I × J. Statistik uji ini diberikan oleh rumus berikut:
\[ G^2 = 2 \sum_i \sum_j n_{ij} \ln\left(\frac{n_{ij}}{\hat{\mu}_{ij}}\right) \]
Di mana: - \(n_{ij}\) adalah frekuensi observasi dalam tabel kontingensi. - \(\hat{\mu}_{ij} = n \cdot p_{i.} \cdot p_{.j}\) adalah frekuensi harapan (ekspektasi), yang dihitung dari total baris dan kolom.
Karakteristik Statistik G²:
- Statistik G² mengikuti distribusi chi-square (χ²)
dengan derajat bebas:
\[ df = (I - 1)(J - 1) \] ##### Perhitungan dengan R
### Frekuensi Harapan:
row_totals <- rowSums(data)
col_totals <- colSums(data)
grand_total <- sum(data)
expected <- outer(row_totals, col_totals) / grand_total
expected## Yes No
## Female 360.7637 84.23629
## Male 260.2363 60.76371### Hitung G² dan p-value:
valid_cells <- data > 0
G2_stat <- 2 * sum(data[valid_cells] *
log(data[valid_cells] / expected[valid_cells]))
df <- (nrow(data) - 1) * (ncol(data) - 1)
p_value <- 1 - pchisq(G2_stat, df)
list(
G2_Statistic = G2_stat,
Degrees_of_Freedom = df,
P_Value = p_value,
Decision = ifelse(p_value < 0.05, "Tolak H0", "Gagal Tolak H0")
)## $G2_Statistic
## [1] 3.628349
##
## $Degrees_of_Freedom
## [1] 1
##
## $P_Value
## [1] 0.05680316
##
## $Decision
## [1] "Gagal Tolak H0"Interpretasi:
Hasil uji chi-kuadrat dengan statistik G² sebesar 3.628349 dan derajat kebebasan 1 menghasilkan nilai p sebesar 0.05680316, yang sedikit lebih besar dari 0,05. Ini berarti bahwa tidak ada perbedaan yang signifikan antara proporsi respons “Ya” dan “Tidak” di antara pria dan wanita dalam data yang dianalisis. Dengan kata lain, berdasarkan hasil uji ini, kita gagal membuktikan adanya perbedaan gender yang signifikan dalam hal proporsi “Ya” dan “Tidak” pada variabel yang diuji.
6.2.3.2 Uji Exact Fisher
Uji Fisher’s Exact digunakan untuk menguji hubungan antara dua variabel kategorikal dalam tabel kontingensi kecil, dimana asumsi Chi-square tidak berlaku karena ukuran sampel yang kecil. Uji ini pertama kali dikembangkan oleh Sir Ronald A. Fisher, seorang ahli statistik, yang menggunakannya dalam penelitian biologi dan medis. Fisher’s Exact Test menjadi sangat penting dalam analisis data yang memiliki frekuensi kecil karena tidak bergantung pada asumsi distribusi normal.
Keunggulan:
- Cocok untuk ukuran sampel kecil.
- Tidak memerlukan asumsi normalitas atau chi-square.
- Memberikan hasil yang lebih akurat dibandingkan uji Chi-square pada data dengan frekuensi kecil.
Keterbatasan: - Perhitungan bisa menjadi sangat berat secara komputasi jika ukuran tabel besar. - Hanya cocok untuk tabel kontingensi kecil (misalnya 2x2 atau 3x3).
Perhitungan dengan R
fisher.test(data)##
## Fisher's Exact Test for Count Data
##
## data: data
## p-value = 0.06165
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.9733491 2.0803837
## sample estimates:
## odds ratio
## 1.423136Interpretasi:
Hasil uji statistik menunjukkan bahwa tidak ada perbedaan signifikan antara pria dan wanita terkait kategori “Ya” dan “Tidak” berdasarkan nilai p chi-square (0,0568) yang sedikit lebih besar dari 0,05. Meskipun ada indikasi perbedaan melalui odds ratio (1,4238) dan nilai p Fisher yang mendekati 0,05 (0,06165), bukti yang ada tidak cukup kuat untuk menolak hipotesis nol.
Kesimpulan
| Uji | Tujuan | Output |
|---|---|---|
| Uji Proporsi | Menguji apakah proporsi kejadian berbeda antara dua kelompok | p-value (signifikansi perbedaan proporsi) |
| Risk Difference (RD) | Mengukur perbedaan absolut risiko antar kelompok | Nilai perbedaan risiko |
| Relative Risk (RR) | Mengukur risiko relatif antar kelompok | Rasio risiko |
| Odds Ratio (OR) | Mengukur peluang kejadian antar kelompok | Rasio odds |
| Uji Chi-Square | Menguji apakah ada hubungan antara dua variabel | p-value (signifikansi hubungan) |
| Uji Fisher’s Exact | Menguji hubungan pada sampel kecil | p-value (signifikansi hubungan) |
6.3 Analisis Residual dalam Tabel Kontingensi
Residual dalam tabel kontingensi digunakan untuk mengidentifkasi sel mana yang menyumbang paling banyak terhadap hubungan antara variabel kategori. Residual mengukur selisih antara frekuensi yang diamati dan frekuensi yang diharapkan berdasarkan model independensi. Jika suatu sel memiliki residual yang besar (positif atau negatif), berarti frekuensi observasi dalam sel tersebut sangat berbeda dari yang diharapkan. Sebaliknya, jika residualnya kecil, berarti nilai observasi mendekati nilai ekspektasi, sehingga sel tersebut tidak banyak menyumbang terhadap hubungan antara variabel (Agresti A., 2013).
Jika residual =0: – Tidak ada perbedaan signifkan antara jumlah observasi yang terjadi di suatu sel dengan jumlah yang diprediksi oleh model independensi. – Artinya, variabel baris dan kolom dalam tabel tidak memiliki hubungan kuat atau tidak menunjukkan pola keterkaitan yang jelas.
Jika residual positif besar: Frekuensi observasi jauh lebih tinggi dari yang diharapkan → menunjukkan hubungan positif yang kuat antara kategori tersebut.
Jika residual negatif besar: – Frekuensi observasi jauh lebih rendah dari yang diharapkan → menunjukkan hubungan negatif atau tidak adanya asosiasi.
Residual kecil (=0) → Tidak ada hubungan antara variabel baris dan kolom.
Residual besar (positif atau negatif) → Ada hubungan yang signifkan antara variabel baris dan kolom.
Jadi, dalam analisis residual tabel kontingensi, nilai residual mendekati nol mengindikasikan bahwa variabel-variabel yang dianalisis tidak memiliki hubungan yang signifkan, dan data sesuai dengan asumsi independensi.
6.3.1 Jenis Residual
- Raw Residual:
\[ r_{ij} = O_{ij} - E_{ij} \]
- Pearson Residual:
\[ r_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}}} \]
- Standardized Residual:
\[ r_{ij}^{*} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}(1 - p_{i.})(1 - p_{.j})}} \]
Perhitungan dengan R
# Data Observasi
data <- matrix(c(20, 10, 30, 20), nrow = 2, byrow = TRUE)
# Hitung nilai ekspektasi
expected <- chisq.test(data)$expected
# Pearson Residual
pearson_residual <- (data - expected) / sqrt(expected)
# Standardized Residual
row_sum <- rowSums(data)
col_sum <- colSums(data)
total_sum <- sum(data)
standardized_residual <- matrix(0, nrow = nrow(data), ncol = ncol(data))
for (i in 1:nrow(data)) {
for (j in 1:ncol(data)) {
standardized_residual[i, j] <- (data[i, j] - expected[i, j]) /
sqrt(expected[i, j] * (1 - row_sum[i] / total_sum) * (1 - col_sum[j] / total_sum))
}
}
# Menampilkan hasil
list(
Pearson_Residual = pearson_residual,
Standardized_Residual = standardized_residual
)## $Pearson_Residual
## [,1] [,2]
## [1,] 0.2886751 -0.3726780
## [2,] -0.2236068 0.2886751
##
## $Standardized_Residual
## [,1] [,2]
## [1,] 0.5962848 -0.5962848
## [2,] -0.5962848 0.59628486.3.2 Deteksi Outlier
Jika nilai residual (khususnya standardized) melebihi ±2, maka sel tersebut dapat dianggap sebagai outlier dalam tabel kontingensi (Agresti, 2018).
Outlier dalam analisis data kategori adalah sel dalam tabel kontingensi yang memiliki nilai residual yang sangat besar, baik positif maupun negatif. Outlier ini menunjukkan bahwa ada kategori yang memiliki frekuensi observasi yang jauh lebih tinggi atau lebih rendah dibandingkan dengan nilai ekspektasi berdasarkan asumsi independensi.
Menggunakan Residual untuk Mendeteksi Outlier
Pearson Residual: \[ e_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}}} \] Jika nilai \(|e_{\{ij\}}| > 2\), maka sel tersebut dianggap sebagai indikasi adanya outlier.
Standardized Residual (Adjusted Residual): \[ r_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}(1 - p_{i\cdot})(1 - p_{\cdot j})}} \] Jika \(|r_{\{ij\}}| > 3\), maka sel tersebut dianggap sebagai outlier signifikan.
Perhitungan dengan R
# Hitung nilai ekspektasi
expected <- chisq.test(data)$expected
# Pearson Residual
pearson_residual <- (data - expected) / sqrt(expected)
# Standardized Residual
row_sum <- rowSums(data)
col_sum <- colSums(data)
total_sum <- sum(data)
standardized_residual <- matrix(0, nrow = nrow(data), ncol = ncol(data))
for (i in 1:nrow(data)) {
for (j in 1:ncol(data)) {
standardized_residual[i, j] <- (data[i, j] - expected[i, j]) /
sqrt(expected[i, j] * (1 - row_sum[i] / total_sum) * (1 - col_sum[j] / total_sum))
}
}
# Menampilkan hasil
list(
Observed = data,
Expected = expected,
Pearson_Residual = pearson_residual,
Standardized_Residual = standardized_residual
)## $Observed
## [,1] [,2]
## [1,] 20 10
## [2,] 30 20
##
## $Expected
## [,1] [,2]
## [1,] 18.75 11.25
## [2,] 31.25 18.75
##
## $Pearson_Residual
## [,1] [,2]
## [1,] 0.2886751 -0.3726780
## [2,] -0.2236068 0.2886751
##
## $Standardized_Residual
## [,1] [,2]
## [1,] 0.5962848 -0.5962848
## [2,] -0.5962848 0.5962848Interpretasi:
Hasil analisis residual menunjukkan bahwa seluruh nilai residual, baik Pearson maupun standardized residual, berada dalam rentang ±1, yang mengindikasikan tidak adanya sel yang menyumbang secara signifikan terhadap ketidaksesuaian antara data observasi dan harapan. Nilai residual terbesar sekitar ±0.60, masih jauh dari ambang ±2 yang umumnya digunakan untuk menandai kontribusi signifikan. Dengan demikian, tidak terdapat penyimpangan mencolok dari ekspektasi, dan distribusi data cenderung sesuai dengan yang diharapkan menurut model independensi.
Bab 7 Tabel Kontingensi Tiga Arah
Analisis tabel kontingensi tiga arah digunakan untuk memahami
hubungan antara tiga variabel kategori secara bersamaan. Sering kali,
hubungan antara dua variabel (misalnya X dan
Y) dapat dipengaruhi oleh variabel ketiga (Z)
yang berperan sebagai variabel pengganggu atau kovariat. Oleh karena
itu, penting membedakan dua jenis tabel:
- Tabel Parsial: Menyajikan hubungan antara
XdanYuntuk setiap level dariZ. Ini memungkinkan analisis hubungan bersyarat dengan mengendalikanZ. - Tabel Marginal: Menyajikan hubungan antara
XdanYdengan mengabaikanZ, yang bisa menyebabkan distorsi interpretasi (misalnya Simpson’s Paradox).
7.1 Tabel Parsial dan Marginal
Tabel marginal diperoleh dengan menjumlahkan frekuensi dari satu variabel. Tabel parsial menyajikan dua variabel untuk setiap level variabel ketiga (Agresti, 2002).
Tabel Kontingensi Tiga Arah
Kita akan membuat tabel kontingensi 3-arah berdasarkan data berikut:
- X (Race): White, Black
- Y (Belief in Afterlife): Yes, Undecided, No
- Z (Gender): Female, Male
# Membuat data utama
data <- data.frame(
Ras = c("White", "White", "Black", "Black"),
Gender = c("Female", "Male", "Female", "Male"),
Yes = c(371, 250, 64, 25),
Undecided = c(49, 45, 9, 5),
No = c(74, 71, 15, 13)
)
# Tampilkan tabel kontingensi utama
knitr::kable(data, caption = "Tabel Kontingensi Tiga Arah: Ras, Gender, dan Kepercayaan terhadap Afterlife")| Ras | Gender | Yes | Undecided | No |
|---|---|---|---|---|
| White | Female | 371 | 49 | 74 |
| White | Male | 250 | 45 | 71 |
| Black | Female | 64 | 9 | 15 |
| Black | Male | 25 | 5 | 13 |
# Tabel Frekuensi Parsial
## Frekuensi Parsial berdasarkan Ras
partial_race <- aggregate(cbind(Yes, Undecided, No) ~ Ras, data = data, sum)
partial_race$Kategori <- paste("Ras:", partial_race$Ras)
## Frekuensi Parsial berdasarkan Gender
partial_gender <- aggregate(cbind(Yes, Undecided, No) ~ Gender, data = data, sum)
partial_gender$Kategori <- paste("Gender:", partial_gender$Gender)
## Frekuensi Parsial berdasarkan Kepercayaan terhadap Afterlife
partial_belief <- data.frame(
Kategori = c("Kepercayaan: Yes", "Kepercayaan: Undecided", "Kepercayaan: No"),
Yes = c(sum(data$Yes), 0, 0),
Undecided = c(0, sum(data$Undecided), 0),
No = c(0, 0, sum(data$No))
)
# Susun ulang kolom agar seragam
partial_race <- partial_race[, c("Kategori", "Yes", "Undecided", "No")]
partial_gender <- partial_gender[, c("Kategori", "Yes", "Undecided", "No")]
# Gabungkan semua menjadi satu tabel
partial_all <- rbind(partial_race, partial_gender, partial_belief)
# Tampilkan hasil akhir frekuensi parsial
knitr::kable(partial_all, caption = "Frekuensi Parsial Gabungan: Ras, Gender, dan Kepercayaan terhadap Afterlife")| Kategori | Yes | Undecided | No |
|---|---|---|---|
| Ras: Black | 89 | 14 | 28 |
| Ras: White | 621 | 94 | 145 |
| Gender: Female | 435 | 58 | 89 |
| Gender: Male | 275 | 50 | 84 |
| Kepercayaan: Yes | 710 | 0 | 0 |
| Kepercayaan: Undecided | 0 | 108 | 0 |
| Kepercayaan: No | 0 | 0 | 173 |
# Tabel Frekuensi Marginal
## Frekuensi marginal baris (total per individu/baris)
data$Total <- rowSums(data[, c("Yes", "Undecided", "No")])
## Frekuensi marginal kolom (total per jawaban)
total_kolom <- colSums(data[, c("Yes", "Undecided", "No", "Total")])
total_kolom <- as.data.frame(t(total_kolom))
rownames(total_kolom) <- "Total"
## Tampilkan tabel frekuensi marginal baris
knitr::kable(data, caption = "Tabel dengan Frekuensi Marginal Baris (Jumlah per Individu)")| Ras | Gender | Yes | Undecided | No | Total |
|---|---|---|---|---|---|
| White | Female | 371 | 49 | 74 | 494 |
| White | Male | 250 | 45 | 71 | 366 |
| Black | Female | 64 | 9 | 15 | 88 |
| Black | Male | 25 | 5 | 13 | 43 |
## Tampilkan tabel frekuensi marginal kolom
knitr::kable(total_kolom, caption = "Frekuensi Marginal Kolom (Jumlah per Kategori Jawaban)")| Yes | Undecided | No | Total | |
|---|---|---|---|---|
| Total | 710 | 108 | 173 | 991 |
7.2 Distribusi Peluang
7.2.1 Peluang Bersama (Joint Probability)
Peluang bersama didefinisikan sebagai: \[ P(Z, X, Y) = \frac{f(Z, X, Y)}{N} \]
7.2.2 Peluang Marginal (Marginal Probability)
Peluang marginal diperoleh dengan menjumlahkan probabilitas bersama:
Peluang terkena kanker \((Y = 1)\): \[ P(Y = 1) = \sum_{Z, X} P(Z, X, Y = 1) \]
Peluang merokok \((X = Y_a)\):
\[ P(X = Y_a) = \sum_{Z, Y} P(Z, X = Y_a, Y) \]
7.3.3 Peluang Bersyarat (Conditional Probability)
Peluang bersyarat didefinisikan sebagai:
\[ P(Z \mid X, Y) = \frac{P(Z, X, Y)}{P(X, Y)} \] ### Perhitungan dengan R
# Tambahkan total masing-masing baris
data$Total <- rowSums(data[, c("Yes", "Undecided", "No")])
# Total keseluruhan
total_all <- sum(data$Total)
# Tampilkan data
knitr::kable(data, caption = "Tabel Data Utama")| Ras | Gender | Yes | Undecided | No | Total |
|---|---|---|---|---|---|
| White | Female | 371 | 49 | 74 | 494 |
| White | Male | 250 | 45 | 71 | 366 |
| Black | Female | 64 | 9 | 15 | 88 |
| Black | Male | 25 | 5 | 13 | 43 |
# Hitung peluang bersama
joint_prob <- data
joint_prob[, c("Yes", "Undecided", "No")] <- round(joint_prob[, c("Yes", "Undecided", "No")] / total_all, 4)
knitr::kable(joint_prob[, c("Ras", "Gender", "Yes", "Undecided", "No")], caption = "Tabel Peluang Bersama")| Ras | Gender | Yes | Undecided | No |
|---|---|---|---|---|
| White | Female | 0.3744 | 0.0494 | 0.0747 |
| White | Male | 0.2523 | 0.0454 | 0.0716 |
| Black | Female | 0.0646 | 0.0091 | 0.0151 |
| Black | Male | 0.0252 | 0.0050 | 0.0131 |
# Peluang marginal untuk setiap kategori jawaban (kolom)
marginal_yes <- sum(data$Yes) / total_all
marginal_undecided <- sum(data$Undecided) / total_all
marginal_no <- sum(data$No) / total_all
# Peluang marginal untuk gender
marginal_gender <- aggregate(Total ~ Gender, data = data, sum)
marginal_gender$Probability <- round(marginal_gender$Total / total_all, 4)
# Peluang marginal untuk ras
marginal_race <- aggregate(Total ~ Ras, data = data, sum)
marginal_race$Probability <- round(marginal_race$Total / total_all, 4)
# Tampilkan
cat("Peluang Marginal: \n")## Peluang Marginal:cat("P(Yes) =", marginal_yes, "\n")## P(Yes) = 0.716448cat("P(Undecided) =", marginal_undecided, "\n")## P(Undecided) = 0.1089808cat("P(No) =", marginal_no, "\n")## P(No) = 0.1745711knitr::kable(marginal_gender[, c("Gender", "Probability")], caption = "Peluang Marginal Berdasarkan Gender")| Gender | Probability |
|---|---|
| Female | 0.5873 |
| Male | 0.4127 |
knitr::kable(marginal_race[, c("Ras", "Probability")], caption = "Peluang Marginal Berdasarkan Ras")| Ras | Probability |
|---|---|
| Black | 0.1322 |
| White | 0.8678 |
# P(Z="Yes" | Ras = "White", Gender = "Female")
numerator <- data[data$Ras == "White" & data$Gender == "Female", "Yes"]
denominator <- data[data$Ras == "White" & data$Gender == "Female", "Total"]
conditional_yes_white_female <- round(numerator / denominator, 4)
cat("P(Yes | White, Female) =", conditional_yes_white_female)## P(Yes | White, Female) = 0.7517.3 Tabel Peluang Bersyarat
Distribusi bersyarat seperti 𝑃(𝑋∣𝑍,), atau 𝑃(𝑌∣𝑋,)), digunakan untuk mengisolasi pengaruh satu variabel terhadap lainnya dalam strata variabel ketiga (Christensen, 1997). \[ P(X = i | Z = k) = \frac{n_{i.k}}{n_{..k}} \] Kesimpulan
- Peluang bersama dihitung dengan membagi frekuensi dengan total populasi.
- Peluang marginal diperoleh dengan menjumlahkan probabilitas bersama.
- Peluang bersyarat dihitung menggunakan rumus Bayes
7.4 Ukuran Asosiasi
7.4.1 Beda Peluang (BP)
\[ BP = P(Y \mid X_1, Z) - P(Y \mid X_2, Z) \]
7.4.2 Risiko Relatif (RR)
\[ RR = \frac{P(Y \mid X_1, Z)}{P(Y \mid X_2, Z)} \]
7.4.3 Odds Ratio (OR)
\[ OR = \frac{ \frac{P(Y \mid X_1, Z)}{1 - P(Y \mid X_1, Z)} }{ \frac{P(Y \mid X_2, Z)}{1 - P(Y \mid X_2, Z)} } \] Ukuran asosiasi dalam tabel kontingensi digunakan untuk mengukur kekuatan hubungan antara dua variabel kategori. Tiga ukuran asosiasi yang umum digunakan adalah Risk Difference (RD), Relative Risk (RR), dan Odds Ratio (OR).
7.4.4 Tabel Kontingensi Parsial
Ukuran Asosiasi
Risk Difference (RD) \[ RD = \frac{a}{a + b} - \frac{c}{c + d} \]
Relative Risk (RR) \[ RR = \frac{\frac{a}{a + b}}{\frac{c}{c + d}} \]
Odds Ratio (OR) \[ OR = \frac{a \cdot d}{b \cdot c} \]
Perhitungan dengan R
kel1 <- data[data$Ras == "White" & data$Gender == "Male", ]
kel2 <- data[data$Ras == "Black" & data$Gender == "Male", ]
p1 <- kel1$Yes / sum(kel1[, c("Yes", "Undecided", "No")])
p2 <- kel2$Yes / sum(kel2[, c("Yes", "Undecided", "No")])
RD <- round(p1 - p2, 4)
RD## [1] 0.1017RR <- round(p1 / p2, 4)
RR## [1] 1.1749odds1 <- p1 / (1 - p1)
odds2 <- p2 / (1 - p2)
OR <- round(odds1 / odds2, 4)
OR## [1] 1.5517Interpretasi:
Individu berjenis kelamin laki-laki dari ras White memiliki kemungkinan 10,17% lebih tinggi untuk percaya pada afterlife dibandingkan dengan laki-laki dari ras Black.
Risiko percaya pada afterlife di antara White Male adalah 17,49% lebih besar dibandingkan Black Male.
Peluang (odds) White Male untuk percaya pada afterlife adalah 55,17% lebih besar dibandingkan dengan Black Male.
7.5 Conditional Independence
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 \mid Z) = P(X \mid Z) \cdot P(Y \mid Z) \]
Dalam bentuk frekuensi:
\[ \frac{n_{ijk}}{n_{k++}} = \frac{n_{i+k}}{n_{k++}} \cdot \frac{n_{+jk}}{n_{k++}} \]
Keterangan: - \(n_{ijk}\): Frekuensi gabungan untuk kategori ke-\(i\), \(j\), dan \(k\) - \(n_{i+k}\): Total frekuensi kategori ke-\(i\) dan \(k\) - \(n_{+jk}\): Total frekuensi kategori ke-\(j\) dan \(k\) - \(n_{k++}\): Total semua frekuensi untuk kategori ke-(k_
Pengujian Conditional Independence
MWnggunakan uji chi-square pada setiap strata.
\[ P(X, Y | Z) \approx P(X | Z)P(Y | Z) \]
7.6 Marginal Y dan X
Distribusi marginal 𝑃(𝑌) atau 𝑃(𝑋) dihitung dengan menjumlahkan nilai probabilitas terhadap variabel lain: \[ P(Y = j) = \sum_{i,k} P(X = i, Y = j, Z = k) \]
7.7 Inferensi Tabel Kontingensi Tiga Arah
Tabel kontingensi tiga arah digunakan untuk menganalisis hubungan antara dua variabel kategorik dengan mempertimbangkan variabel kontrol.
Tabel ini terdiri dari beberapa tabel parsial (2x2) untuk setiap tingkat \(Z\), serta tabel marginal yang mengabaikan \(Z\). Ukuran asosiasi yang digunakan adalah odds ratio:
\[ OR = \frac{a \cdot d}{b \cdot c} \]
Jika odds ratio parsial relatif konstan, kita dapat menghitung odds ratio bersama menggunakan estimasi Mantel-Haenszel.
7.7.1 Independensi Bersyarat
Independensi bersyarat adalah konsep penting dalam analisis tabel kontingensi tiga arah. Ini merujuk pada kondisi di mana dua variabel, X dan Y independen dalam setiap level variabel ketiga, Z. Pengujian independensi bersyarat dilakukan dengan metode statistik seperti uji Cochran-Mantel-Haenszel (CMH).
Independensi Bersyarat: Dua variabel, X dan Y dikatakan independen bersyarat terhadap variabel ketiga, Z, jika rasio odds mereka dalam setiap strata Z sama dengan 1.
\[ OR(X, Y | Z) = 1 \]
Artinya, setelah mengendalikan pengaruh \(Z\), tidak ada hubungan antara \(X\) dan \(Y\) dalam setiap strata.
Hal yang Perlu Diperhatikan
- Meskipun \(X\) dan \(Y\) independen dalam setiap strata \(Z\), hal ini tidak berarti bahwa mereka juga independen dalam tabel marginal (keseluruhan data tanpa mempertimbangkan \(Z\)). Fenomena ini dikenal sebagai Paradox Simpson.
- Asumsi independensi bersyarat umum digunakan dalam studi epidemiologi, ilmu sosial, dan model inferensial lainnya yang menggunakan tabel kontingensi struktural.
7.7.2 Uji Statistik Independensi Bersyarat
Metode Cochran-Mantel-Haenszel
Tujuan Uji CMH digunakan untuk menguji hubungan antara dua variabel kategori dengan mempertimbangkan efek dari variabel perancu (confounder). Uji ini berguna untuk:
- Menguji hubungan yang dikendalikan oleh faktor perancu dalam tabel kontingensi berlapis (stratified contingency tables).
- Mengendalikan variabel perancu antara dua variabel kategori dengan mempertimbangkan efek dari variabel ketiga.
- Mengatasi bias akibat confounding, misalnya dalam studi epidemiologi dan eksperimen sosial.
Ide Dasar Uji CMH berangkat dari konsep tabel kontingensi berlapis (stratified 2 × 2 tables), yang menggabungkan informasi dari strata berbeda untuk menguji efek bersama variabel perancu.
Sebagai contoh, jika kita ingin menguji hubungan antara merokok (\(X\)) dan kanker paru-paru (\(Y\)), namun ingin mengontrol efek dari usia (\(Z\)), maka kita membentuk tabel 2×2 untuk setiap strata usia. CMH kemudian menguji hubungan antara \(X\) dan \(Y\) setelah mengontrol efek \(Z\).
Hipotesis
- Hipotesis nol (\(H_0\)): \(\theta_{XY|k} = 1\) untuk setiap \(k = 1, 2, ..., K\)
- Hipotesis alternatif (\(H_1\)): \(\theta_{XY|k} \ne 1\) untuk paling sedikit satu \(k\)
Statistik uji Cochran-Mantel-Haenszel (CMH) dirumuskan sebagai:
\[ CMH = \frac{\left[ \sum_k (n_{11k} - \mu_{11k}) \right]^2}{\sum_k \text{var}(n_{11k})} \]
- \(n_{11k}\): nilai frekuensi sel baris 1 kolom 1 pada tabel parsial ke-\(k\)
- \(\mu_{11k}\): nilai ekspektasi sel baris 1 kolom 1 pada tabel parsial ke-\(k\), dihitung dengan rumus:
\[ \mu_{11k} = E(n_{11k}) = \frac{n_{1 \cdot k} \cdot n_{\cdot 1 k}}{n_{\cdot \cdot k}} \] - Varian dari \(n_{11k}\) diberikan oleh: \[ \text{var}(n_{11k}) = \frac{n_{1 \cdot k} \cdot n_{2 \cdot k} \cdot n_{\cdot 1 k} \cdot n_{\cdot 2 k}}{n_{\cdot \cdot k}^2 (n_{\cdot \cdot k} - 1)} \] - Statistik CMH mengikuti distribusi Chi-square dengan derajat kebebasan (\(df\)) = 1. - Tolak \(H_0\) jika \(CMH > \chi^2_{(1)}\) atau p-value < α.
Perhitungan dengan R
data_cmh <- data.frame(
Ras = rep(c("White", "White", "Black", "Black"), each = 1),
Gender = c("Female", "Male", "Female", "Male"),
Yes = c(371, 250, 64, 25),
No = c(74, 71, 15, 13)
)
# Membuat array 2x2x2 (Gender x Belief x Ras)
cmh_array <- array(c(
371, 74, # White Female: Yes, No
250, 71, # White Male: Yes, No
64, 15, # Black Female: Yes, No
25, 13 # Black Male: Yes, No
), dim = c(2, 2, 2),
dimnames = list(
Gender = c("Female", "Male"),
Belief = c("Yes", "No"),
Ras = c("White", "Black")
))
cmh_array## , , Ras = White
##
## Belief
## Gender Yes No
## Female 371 250
## Male 74 71
##
## , , Ras = Black
##
## Belief
## Gender Yes No
## Female 64 25
## Male 15 13# Melakukan uji CMH
mantelhaen.test(cmh_array)##
## Mantel-Haenszel chi-squared test with continuity correction
##
## data: cmh_array
## Mantel-Haenszel X-squared = 5.5778, df = 1, p-value = 0.01819
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 1.085101 2.120697
## sample estimates:
## common odds ratio
## 1.516961Interpretasi:
Berdasarkan uji Mantel-Haenszel terhadap data yang menggabungkan sikap pemilih berdasarkan ras (White vs Black) dan gender (Female vs Male), ditemukan bahwa terdapat asosiasi yang signifikan antara ras dan pilihan suara (Yes, Undecided, No) setelah mengontrol variabel gender (p-value = 0.018). Hasil ini menunjukkan bahwa secara keseluruhan, ras memengaruhi pilihan suara bahkan ketika perbedaan gender diperhitungkan. Nilai common odds ratio sebesar 1.52 mengindikasikan bahwa individu dari satu ras (kemungkinan besar White, jika diurutkan sebagai referensi pertama) memiliki sekitar 1,5 kali lebih besar peluang untuk memilih “Yes” dibandingkan individu dari ras lainnya (Black), setelah memperhitungkan gender. Confidence interval yang tidak mencakup 1 (1.09–2.12) memperkuat kesimpulan bahwa asosiasi ini statistik signifikan.
7.7.3 Odds Ratio Bersama
Penaksir (Khusus Tabel \(2 \times 2 \times K\))
- Dalam tabel kontingensi \(2 \times 2 \times K\), terdapat \(K\) tabel parsial, sehingga ada sebanyak \(K\) odds ratio bersyarat.
- Jika nilai-nilai odds ratio bersyarat relatif sama (tidak berbeda secara ekstrem) dan memiliki arah yang sama, maka kita dapat menentukan sebuah nilai tunggal untuk odds ratio yang disebut odds ratio bersama.
- Odds ratio bersama ditaksir oleh statistik Mantel-Haenszel.
Rumus Odds Ratio Bersama
Odds ratio bersama ditaksir menggunakan rumus: \[ \hat{\theta}_{MH} = \frac{\sum_{k=1}^{K} \frac{n_{11k}n_{22k}}{n_{++k}}}{\sum_{k=1}^{K} \frac{n_{12k}n_{21k}}{n_{++k}}} \]
Keterangan: - \(n_{11k}\): Frekuensi sel baris 1 kolom 1 pada tabel parsial ke-\(k\) - \(n_{12k}\): Frekuensi sel baris 1 kolom 2 pada tabel parsial ke-\(k\) - \(n_{21k}\): Frekuensi sel baris 2 kolom 1 pada tabel parsial ke-\(k\) - \(n_{22k}\): Frekuensi sel baris 2 kolom 2 pada tabel parsial ke-\(k\) - \(n_{++k}\): Total observasi dalam tabel parsial ke-\(k\)
Standard Error untuk log Odds Ratio Bersama
Standard error untuk log odds ratio bersama dihitung dengan rumus:
\[ \hat{\sigma}^2[\log(\hat{\theta}_{MH})] = \frac{\sum (n_{11k} + n_{12k})(n_{11k}n_{22k}) / n_{..k}^2}{2 \left( \sum n_{11k}n_{22k} / n_{..k} \right)^2} + \frac{\sum \left[(n_{11k} + n_{22k})(n_{11k} + n_{12k}) + (n_{12k} + n_{21k})(n_{11k} + n_{22k})\right] / n_{..k}^2}{ 2 \left( \sum n_{11k}n_{22k} / n_{..k} \right) \left( \sum n_{12k}n_{21k} / n_{..k} \right) } + \frac{\sum (n_{12k} + n_{21k})(n_{12k}n_{21k}) / n_{..k}^2}{2 \left( \sum n_{12k}n_{21k} / n_{..k} \right)^2} \]
Interval Kepercayaan log Odds Ratio Bersama
\[ \log(\hat{\theta}_{MH}) \pm Z_{\alpha/2} \hat{\sigma}[\log(\hat{\theta}_{MH})] \] ### Perhitungan dengan R
# Uji CMH lagi untuk mendapatkan log odds ratio dan CI
cmh_result <- mantelhaen.test(cmh_array)
# Log odds ratio (diambil dari statistik)
log_or <- log(cmh_result$estimate)
log_or## common odds ratio
## 0.416709# Standard Error (SE) dari log odds ratio bisa dihitung dari CI:
# log(CI upper) - log(CI lower) = 2 * 1.96 * SE
ci_log <- log(cmh_result$conf.int)
se_log_or <- (ci_log[2] - ci_log[1]) / (2 * 1.96)
se_log_or## [1] 0.1709367# Interval Kepercayaan 95% dari log odds ratio
lower_log_or <- log_or - 1.96 * se_log_or
upper_log_or <- log_or + 1.96 * se_log_or
# Interval Kepercayaan 95% dari odds ratio (eksp dari log CI)
ci_or <- exp(c(lower_log_or, upper_log_or))
# Tampilkan hasil
list(
"Log Odds Ratio" = log_or,
"Standard Error (log OR)" = se_log_or,
"95% CI (log OR)" = c(lower_log_or, upper_log_or),
"Odds Ratio" = cmh_result$estimate,
"95% CI (OR)" = ci_or
)## $`Log Odds Ratio`
## common odds ratio
## 0.416709
##
## $`Standard Error (log OR)`
## [1] 0.1709367
##
## $`95% CI (log OR)`
## common odds ratio common odds ratio
## 0.08167294 0.75174497
##
## $`Odds Ratio`
## common odds ratio
## 1.516961
##
## $`95% CI (OR)`
## common odds ratio common odds ratio
## 1.085101 2.120697Interpretasi:
Hasil uji Mantel-Haenszel menunjukkan bahwa terdapat asosiasi yang signifikan antara ras dan preferensi suara setelah mengendalikan gender, dengan odds ratio sebesar 1.52 dan interval kepercayaan 95% antara 1.09 hingga 2.12, yang menunjukkan bahwa individu dari ras tertentu (kemungkinan White) memiliki peluang lebih besar untuk memilih “Yes” dibandingkan ras lainnya (Black). Dalam bentuk logaritmik, log odds ratio sebesar 0.42 dengan standard error 0.17, dan interval kepercayaan 95% untuk log odds ratio berada antara 0.082 hingga 0.752. Karena interval ini tidak mencakup nol, ini semakin memperkuat bukti bahwa asosiasi tersebut signifikan secara statistik. Dengan kata lain, bahkan setelah memperhitungkan pengaruh gender, ras tetap menjadi faktor yang berpengaruh terhadap pilihan suara.
7.7.4 Uji Homogenitas (Breslow-Day)
Homogenitas Asosiasi dalam Tabel Kontingensi Tiga Arah
Definisi Asosiasi Homogen
- Asosiasi homogen terjadi jika odds ratio pada setiap tabel parsial bernilai sama:
\[ \theta_{XY(1)} = \theta_{XY(2)} = \cdots = \theta_{XY(K)} \]
Jika odds ratio konstan di semua strata variabel kontrol (\(Z\)), maka tidak ada interaksi antara variabel \(X\) dan \(Y\) dalam pengaruhnya terhadap \(Z\).
Jika odds ratio berbeda-beda antar level \(Z\), maka terdapat interaksi antara \(X\) dan \(Y\) dalam hubungannya dengan \(Z\).
Pengujian Homogenitas dengan Statistik Breslow-Day
Hipotesis
Hipotesis nol (\(H_0\)) : \(\theta_{XY(1)} = \theta_{XY(2)} = \cdots = \theta_{XY(K)}\)
(odds ratio sama di seluruh strata)Hipotesis alternatif (\(H_1\)) : Setidaknya satu odds ratio yang berbeda
Statistik Uji Breslow-Day
Statistik uji Breslow-Day (BD) digunakan untuk menguji homogenitas odds ratio:
\[ X^2_{BD} = \sum_{k=1}^{K} \frac{(a_j - \hat{a}_j)^2}{\text{Var}(\hat{a}_j | \hat{\theta}_{OR, MH})} \]
Dengan: - \(a_j\) = adalah jumlah kasus terpapar yang diamati dalam strata \(j\) - \(\hat{a}_j\) = adalah jumlah kasus terpapar yang diharapkan berdasarkan hipotesis nol - \(\text{Var}(\hat{a}_j | \hat{\theta}_{OR, MH})\) adalah varians dari \(\hat{a}_j\)
- Statistik BD mengikuti distribusi Chi-square dengan derajat kebebasan \(df = K - 1\).
- Jika nilai BD lebih besar dari nilai kritis \(\chi^2_{(K-1)}\), maka tolak \(H_0\). Jika tidak, gagal menolak \(H_0\).
Perhitungan dengan R
# Data input
data_white <- matrix(c(371, 74, 250, 71), nrow = 2, byrow = TRUE,
dimnames = list(Gender = c("Female", "Male"),
Belief = c("Yes", "No")))
data_black <- matrix(c(64, 15, 25, 13), nrow = 2, byrow = TRUE,
dimnames = list(Gender = c("Female", "Male"),
Belief = c("Yes", "No")))
data_white## Belief
## Gender Yes No
## Female 371 74
## Male 250 71data_black## Belief
## Gender Yes No
## Female 64 15
## Male 25 13# Ekstraksi elemen dari masing-masing tabel
a1 <- 371; b1 <- 74; c1 <- 250; d1 <- 71; n1 <- a1 + b1 + c1 + d1
a2 <- 64; b2 <- 15; c2 <- 25; d2 <- 13; n2 <- a2 + b2 + c2 + d2
# Odds Ratio per strata
or_white <- (a1 * d1) / (b1 * c1)
or_black <- (a2 * d2) / (b2 * c2)
or_white## [1] 1.423838or_black## [1] 2.218667# Komponen CMH
num_cmh <- (a1 * d1 / n1) + (a2 * d2 / n2)
den_cmh <- (b1 * c1 / n1) + (b2 * c2 / n2)
or_cmh <- num_cmh / den_cmh
or_cmh## [1] 1.516961# E_ij (expected counts), per strata berdasarkan OR_CMH
# Kita gunakan formula perkiraan expected count dari OR gabungan
# Fungsi untuk menghitung expected a_ij untuk tiap strata
expected_a <- function(a, b, c, d, or_cmh) {
n <- a + b + c + d
row1 <- a + b
col1 <- a + c
A_hat <- (row1 * col1 * or_cmh) / (col1 * or_cmh + (n - col1))
return(A_hat)
}
# Expected cell a untuk masing-masing strata
e1 <- expected_a(a1, b1, c1, d1, or_cmh)
e2 <- expected_a(a2, b2, c2, d2, or_cmh)
# Varians approx untuk tiap strata
var1 <- 1 / e1 + 1 / (a1 + b1 - e1) + 1 / (a1 + c1 - e1) + 1 / (n1 - a1 - b1 - c1 + e1)
var2 <- 1 / e2 + 1 / (a2 + b2 - e2) + 1 / (a2 + c2 - e2) + 1 / (n2 - a2 - b2 - c2 + e2)
# Statistik uji Breslow-Day
bd_stat <- ((a1 - e1)^2 / var1) + ((a2 - e2)^2 / var2)
bd_stat## [1] 8297.652# Derajat bebas = jumlah strata - 1 = 1
p_value <- 1 - pchisq(bd_stat, df = 1)
p_value## [1] 0Interpretasi:
Hasil uji Breslow-Day menunjukkan nilai statistik sebesar 8297.65 dengan p-value = 0, yang jauh di bawah batas signifikansi umum (misalnya 0.05). Ini menunjukkan bahwa terdapat ketidakhomogenan odds ratio antar strata (dalam hal ini, antar gender). Artinya, efek ras terhadap preferensi suara tidak konsisten antara laki-laki dan perempuan—misalnya, pengaruh ras terhadap pilihan suara bisa lebih kuat pada satu gender dibanding yang lain. Oleh karena itu, asumsi homogenitas yang menjadi dasar uji Mantel-Haenszel tidak terpenuhi, dan hasil CMH sebaiknya ditindaklanjuti dengan analisis yang mempertimbangkan efek interaksi atau dilakukan secara terpisah per strata.
Bab 8 Generalized Linear Model (GLM)
Generalized Linear Model (GLM) merupakan pengembangan dari regresi linear klasik. GLM memungkinkan kita menganalisis data dengan distribusi respons yang tidak harus normal dan memungkinkan adanya hubungan non-linear antara prediktor dan rata-rata respons.
GLM memiliki tiga komponen utama:
- Distribusi dari exponential family yang digunakan untuk variabel respons.
- Fungsi link, yaitu fungsi yang menghubungkan rata-rata respons dengan kombinasi linear dari prediktor.
- Fungsi linear prediktor: \(\eta = X \beta\)
8.1 Exponential Family
Suatu distribusi termasuk dalam exponential family jika dapat ditulis dalam bentuk:
\[ f(y; \theta, \phi) = \exp\left\{ \frac{y \theta - b(\theta)}{a(\phi)} + c(y, \phi) \right\} \]
Beberapa distribusi yang termasuk dalam exponential family antara lain:
- Distribusi Normal
- Distribusi Binomial
- Distribusi Poisson
- Distribusi Gamma
Contoh: Distribusi Binomial
Fungsi peluang distribusi binomial:
\[ P(Y = y) = \binom{n}{y} \pi^y (1 - \pi)^{n - y} \]
Dapat ditulis kembali dalam bentuk exponential family:
\[ P(Y = y) = \exp\left\{ \log \binom{n}{y} + y \log\left( \frac{\pi}{1 - \pi} \right) + n \log(1 - \pi) \right\} \]
Dengan substitusi: - \(\theta = \log\left(\frac{\pi}{1 - \pi}\right)\) - \(b(\theta) = -n \log(1 - \pi)\) - \(a(\phi) = 1\)
Maka, distribusi binomial termasuk dalam exponential family.
8.2 Model Regresi Logistik
Model regresi logistik serupa dengan regresi linear, namun hasil prediksinya dibatasi dalam rentang nilai tertentu (biasanya 0 hingga 1). Model ini digunakan untuk memodelkan variabel respons biner dengan pendekatan probabilistik menggunakan fungsi aktivasi sigmoid.
Ciri-ciri model regresi logistik: - Prediksi berada dalam rentang 0 hingga 1. - Kurva yang dihasilkan berbentuk S (sigmoid).
Aplikasi umum model logistik: - Menentukan probabilitas penyakit jantung berdasarkan faktor risiko. - Menyaring pelamar untuk beasiswa atau program pendidikan berdasarkan data nilai dan latar belakang. - Menentukan apakah email merupakan spam atau tidak berdasarkan teks email.
Kelebihan Model Regresi Logistik
- Mudah dipahami dan diinterpretasikan: Model ini cocok untuk data kategorik.
- Menghindari prediksi di luar batas logis: Tidak seperti regresi linear, regresi logistik menjaga prediksi tetap dalam [0,1].
- Persyaratan asumsi yang lebih fleksibel: Tidak membutuhkan hubungan linear antara prediktor dan output.
Fungsi Sigmoid
Fungsi sigmoid didefinisikan sebagai:
\[ f(z) = \frac{1}{1 + e^{-z}} \]
Interpretasi hasil fungsi sigmoid:
- Jika \(f(z) > 0.5\), maka output diklasifikasikan sebagai 1 (positif).
- Jika \(f(z) < 0.5\), maka output diklasifikasikan sebagai 0 (negatif).
Sebagai contoh, jika hasil sigmoid adalah 0.65, maka ada peluang sebesar 65% bahwa peristiwa yang dimodelkan akan terjadi.
Perhitungan dengan R
library(tidyverse)## Warning: package 'tidyverse' was built under R version 4.3.3## Warning: package 'ggplot2' was built under R version 4.3.3## Warning: package 'tibble' was built under R version 4.3.3## Warning: package 'tidyr' was built under R version 4.3.3## Warning: package 'readr' was built under R version 4.3.3## Warning: package 'purrr' was built under R version 4.3.3## Warning: package 'dplyr' was built under R version 4.3.3## Warning: package 'stringr' was built under R version 4.3.3## Warning: package 'forcats' was built under R version 4.3.3## Warning: package 'lubridate' was built under R version 4.3.3## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.0 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(broom)## Warning: package 'broom' was built under R version 4.3.3# Buat data mentah dalam format data.frame
# Data berdasarkan tabel
data <- data.frame(
Race = c("White", "White", "White", "White",
"Black", "Black", "Black", "Black"),
Gender = c("Female", "Male", "Female", "Male",
"Female", "Male", "Female", "Male"),
Belief = c("Yes", "Yes", "Undecided", "Undecided",
"Yes", "Yes", "Undecided", "Undecided"),
Count = c(371, 250, 49, 45, 64, 25, 9, 5)
)
# Tambah data "No"
data_no <- data.frame(
Race = c("White", "White", "Black", "Black"),
Gender = c("Female", "Male", "Female", "Male"),
Belief = "No",
Count = c(74, 71, 15, 13)
)
# Gabungkan semua
data_full <- bind_rows(data, data_no)
# Ubah tipe
data_full <- data_full %>%
mutate(Race = factor(Race),
Gender = factor(Gender),
Belief = factor(Belief, levels = c("No", "Undecided", "Yes")))
data_full# Buat data biner (Yes = 1, selain itu = 0)
data_binary <- data_full %>%
mutate(Yes = ifelse(Belief == "Yes", 1, 0)) %>%
group_by(Race, Gender, Yes) %>%
summarise(Count = sum(Count), .groups = "drop")
data_binary# Model Regresi logistik : Yes ~ Race + Gender
model <- glm(cbind(Count, ave(Count, Race, Gender, FUN = sum) - Count) ~ Race + Gender + Yes,
data = data_binary, family = binomial)
summary(model)##
## Call:
## glm(formula = cbind(Count, ave(Count, Race, Gender, FUN = sum) -
## Count) ~ Race + Gender + Yes, family = binomial, data = data_binary)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -9.269e-01 1.496e-01 -6.196 5.8e-10 ***
## RaceWhite 3.631e-16 1.475e-01 0.000 1
## GenderMale -1.934e-16 1.015e-01 0.000 1
## Yes 1.854e+00 9.967e-02 18.599 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 401.148 on 7 degrees of freedom
## Residual deviance: 17.154 on 4 degrees of freedom
## AIC: 68.114
##
## Number of Fisher Scoring iterations: 3# Prediksi Probabilitas
data_pred <- expand.grid(
Race = c("White", "Black"),
Gender = c("Female", "Male"),
Yes = 1
)
data_pred$prob <- predict(model, newdata = data_pred, type = "response")
ggplot(data_pred, aes(x = interaction(Race, Gender), y = prob, fill = Race)) +
geom_col(position = "dodge") +
labs(title = "Probabilitas Percaya Kehidupan Setelah Mati",
x = "Kelompok", y = "Probabilitas") +
ylim(0,1)
# Kurva Sigmoid
sigmoid <- function(x) {
1 / (1 + exp(-x))
}
df_curve <- tibble(x = seq(-6, 6, length.out = 100),
prob = sigmoid(x))
ggplot(df_curve, aes(x = x, y = prob)) +
geom_line(color = "blue", size = 1.2) +
labs(title = "Kurva Sigmoid", x = "Nilai logit", y = "Probabilitas")## 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.
Interpretasi:
Model regresi logistik digunakan untuk memprediksi peluang seseorang percaya pada afterlife berdasarkan ras dan gender. Koefisien regresi (dengan exp(coef)) menunjukkan odds ratio relatif terhadap referensi. Kelompok White Female memiliki probabilitas tertinggi percaya. Kelompok Black Male memiliki peluang paling kecil percaya. Kurva sigmoid menggambarkan transisi probabilitas dari 0 ke 1 saat nilai logit meningkat.
Hasil regresi logistik biner menunjukkan bahwa variabel Yes (indikator untuk memilih “Yes”) memiliki pengaruh signifikan terhadap peluang memilih, dengan koefisien 1.854 dan p-value < 2e-16, yang berarti semakin tinggi nilai pada variabel ini, semakin besar kemungkinan responden memilih “Yes”. Sebaliknya, variabel Race (White) dan Gender (Male) tidak menunjukkan pengaruh signifikan terhadap peluang memilih “Yes”, karena keduanya memiliki nilai p sebesar 1. Ini menunjukkan bahwa setelah mengendalikan variabel lainnya, ras dan gender tidak memiliki pengaruh signifikan secara statistik, dan hanya preferensi terhadap pilihan “Yes” yang benar-benar memengaruhi model. Model ini juga memiliki peningkatan kecocokan yang signifikan, terlihat dari penurunan deviance dari 401.15 menjadi 17.15, menunjukkan model cukup baik dalam menjelaskan variasi data.
8.3 Model Regresi Poisson
Model regresi Poisson digunakan ketika variabel respons berbentuk data cacah (count data), yaitu bilangan bulat non-negatif. Model ini termasuk dalam Generalized Linear Model (GLM) dengan asumsi bahwa variabel respons mengikuti distribusi Poisson.
Distribusi Poisson
Distribusi Poisson memiliki fungsi probabilitas sebagai berikut:
\[ P(Y = y) = \frac{e^{-\lambda} \lambda^y}{y!} \]
Distribusi ini dapat ditulis ke dalam bentuk exponential family:
\[ f(y; \theta) = \exp \left\{ y \log(\lambda) - \lambda - \log(y!) \right\} \]
Dengan substitusi: - \(\theta = \log(\lambda)\) - \(b(\theta) = e^\theta = \lambda\) - \(\phi = 1\) - \(c(y, \phi) = -\log(y!)\)
Sehingga, distribusi Poisson termasuk dalam keluarga exponential.
Fungsi Link
Fungsi link kanonik yang digunakan untuk distribusi Poisson adalah fungsi logaritma:
\[ g(\mu) = \log(\mu) \]
Modelnya dapat dituliskan sebagai:
\[ \log(\mu_i) = \mathbf{x}_i^T \beta \]
dan fungsi link inversenya adalah:
\[ \mu_i = \exp(\mathbf{x}_i^T \beta) \]
Estimasi Parameter
Parameter \(\beta\) dalam model regresi Poisson dapat diestimasi menggunakan metode Maximum Likelihood Estimation (MLE).
Fungsi log-likelihood untuk model ini adalah:
\[ l(\beta) = \sum_{i=1}^n \left[ y_i \mathbf{x}_i^T \beta - \exp(\mathbf{x}_i^T \beta) - \log(y_i!) \right] \]
Estimasi nilai \(\beta\) dapat diperoleh melalui pendekatan numerik seperti iterasi metode Newton-Raphson.
Perhitungan dengan R
# Ubah ke faktor
data_full <- data_full %>%
mutate(
Race = factor(Race),
Gender = factor(Gender),
Belief = factor(Belief, levels = c("No", "Undecided", "Yes"))
)
data_full# Model Poisson: memodelkan jumlah orang berdasarkan faktor
poisson_model <- glm(Count ~ Race + Gender + Belief, family = poisson(link = "log"), data = data_full)
summary(poisson_model)##
## Call:
## glm(formula = Count ~ Race + Gender + Belief, family = poisson(link = "log"),
## data = data_full)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 2.59753 0.11452 22.683 < 2e-16 ***
## RaceWhite 1.88174 0.09379 20.063 < 2e-16 ***
## GenderMale -0.35276 0.06452 -5.467 4.57e-08 ***
## BeliefUndecided -0.47116 0.12264 -3.842 0.000122 ***
## BeliefYes 1.41197 0.08479 16.653 < 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: 1265.541 on 11 degrees of freedom
## Residual deviance: 14.141 on 7 degrees of freedom
## AIC: 90.454
##
## Number of Fisher Scoring iterations: 4# Tambahkan prediksi ke dalam data
data_full <- data_full %>%
mutate(
Predicted = predict(poisson_model, type = "response")
)
# Visualisasi
ggplot(data_full, aes(x = interaction(Race, Gender, Belief), y = Count)) +
geom_col(fill = "steelblue", alpha = 0.6) +
geom_point(aes(y = Predicted), color = "red", size = 3) +
theme_minimal() +
labs(title = "Observed vs Predicted Counts (Regresi Poisson)",
x = "Group (Race-Gender-Belief)",
y = "Count")
# Simulasi residual menggunakan DHARMa
library(DHARMa)## Warning: package 'DHARMa' was built under R version 4.3.3## Warning in check_dep_version(): ABI version mismatch:
## lme4 was built with Matrix ABI version 1
## Current Matrix ABI version is 0
## Please re-install lme4 from source or restore original 'Matrix' package## This is DHARMa 0.4.7. For overview type '?DHARMa'. For recent changes, type news(package = 'DHARMa')sim <- simulateResiduals(poisson_model)
plot(sim)
# Uji Overdispersion secara manual
resid_dev <- sum(residuals(poisson_model, type = "deviance")^2)
df_resid <- df.residual(poisson_model)
dispersion_stat <- resid_dev / df_resid
dispersion_stat## [1] 2.020133Interpretasi:
Model regresi Poisson menunjukkan bahwa semua variabel prediktor—ras (White), gender (Male), dan keyakinan (Belief)—berpengaruh signifikan terhadap jumlah responden (Count), dengan semua p-value < 0.001. Secara spesifik, individu berkulit putih dan yang memiliki keyakinan “Yes” cenderung memiliki jumlah yang lebih tinggi, sedangkan laki-laki dan responden yang “Undecided” cenderung memiliki jumlah lebih rendah. Model ini memiliki kecocokan yang sangat baik terhadap data, terlihat dari penurunan deviance dari 1265.54 menjadi 14.14, serta nilai AIC yang relatif rendah (90.45). Namun, hasil uji overdispersi menunjukkan nilai statistik dispersi sebesar 2.02, yang berarti model mengalami overdispersion, yaitu varians residual lebih besar dari yang diharapkan dalam model Poisson ideal (nilai ideal = 1). Hal ini disarankan untuk ditangani dengan model alternatif seperti negative binomial regression. Visualisasi juga menunjukkan bahwa prediksi model cukup sesuai dengan data observasi, meskipun adanya overdispersion menandakan kemungkinan ada pola yang tidak sepenuhnya ditangkap oleh model saat ini.
Bab IX Inferensi GLM
Dalam Generalized Linear Model (GLM), penting untuk memahami bagaimana ekspektasi dan varians dari estimator bekerja, karena ini merupakan dasar dari inferensi statistik dan pengujian hipotesis seperti Wald test, Likelihood Ratio Test, serta pembentukan interval kepercayaan.
1. Ekspektasi Estimator
Ekspektasi menjelaskan apakah suatu estimator bersifat tidak bias. Estimator dikatakan tidak bias jika:
\[ E[\hat{\beta}] = \beta \]
Pada GLM, estimasi Maximum Likelihood Estimator (MLE) dari \(\hat{\beta}\) bersifat asymptotically unbiased.
2. Varians Estimator
Varians menggambarkan seberapa presisi sebuah estimator. Varians dari estimator \(\hat{\beta}\) dapat dihampiri dengan:
\[ \text{Var}(\hat{\beta}) \approx (X^T W X)^{-1} \]
dengan \(W\) adalah matriks bobot yang tergantung pada distribusi dan fungsi link.
Distribusi Asimptotik Estimator
Untuk ukuran sampel besar:
\[ \hat{\beta} \sim \mathcal{N}(\beta, \text{Var}(\hat{\beta})) \]
Distribusi ini digunakan untuk:
- Uji Wald
- Confidence Interval
- P-value
Varians dalam GLM Tidak Konstan
Tidak seperti regresi linear biasa (OLS) yang mengasumsikan varian konstan:
\[ \text{Var}(Y_i) = \sigma^2 \]
Dalam GLM, varian model dinyatakan sebagai:
\[ \text{Var}(Y_i) = \phi V(\mu_i) \]
dengan: - \(\phi\) = parameter dispersi - \(V(\mu)\) = fungsi varian
Contoh: Poisson: \(V(\mu) = \mu\) dan Binomial: \(V(\mu) = \mu(1 - \mu)\)
Kesimpulan
- Ekspektasi digunakan untuk mengetahui bias dari suatu estimasi.
- Varians digunakan untuk mengukur presisi dan membangun uji statistik.
- Estimator MLE bersifat asimptotik dan mengikuti distribusi normal untuk sampel besar.
- Varians dalam GLM bergantung pada bentuk fungsi distribusi eksponensial dari data.
9.1 Ekspektasi dan Varians
9.1.1 Ekspektasi
Dari fungsi momen:
\[ E(Y) = \int y f(y; \theta) dy = \mu \]
Untuk keluarga eksponensial:
\[ \log f(y; \theta) = a(y) + b(\theta)y + c(\theta) \]
Atau:
\[ \log f(y; \theta) = y\theta - b(\theta) + c(y) \]
Turunan pertama dari fungsi log-likelihood (score function):
\[ U(\theta) = \frac{\partial \ell}{\partial \theta} = y - b'(\theta) \]
Dengan ekspektasi:
\[ E[U(\theta)] = \mu - b'(\theta) = 0 \Rightarrow \mu = b'(\theta) \]
9.1.2 Varians
Turunan kedua:
\[ \frac{\partial^2 \ell}{\partial \theta^2} = -b''(\theta) \]
Sehingga:
\[ \text{Var}(Y) = b''(\theta) = \phi V(\mu) \]
9.2 Penaksiran Parameter (Maximum Likelihood)
Maximum Likelihood Estimation (MLE)
MLE dilakukan dengan memaksimalkan fungsi likelihood. Langkah dasar:
- Turunan pertama = 0
- Turunan kedua < 0
Karena GLM tidak selalu memiliki bentuk eksplisit, pendekatan numerik digunakan.
Metode Newton-Raphson
Menggunakan: - Score vector (gradien) - Hessian matrix
Iterasi:
\[ \beta^{(t+1)} = \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \]
Fisher Scoring
Modifikasi Newton-Raphson dengan mengganti Hessian menjadi matriks informasi Fisher.
IRLS (Iteratively Reweighted Least Square)
Modifikasi Fisher Scoring yang hasilnya menyerupai metode Least Square.
Implementasi Newton-Raphson
- Statistik score ke-\(j\):
\[ U_j(\beta) = \frac{\partial \log L(\beta)}{\partial \beta_j} \]
- Turunan kedua (Hessian):
\[ H_{jk}(\beta) = \frac{\partial^2 \log L(\beta)}{\partial \beta_j \partial \beta_k} \]
- Taylor expansion:
\[ U(\beta^*) \approx U(\beta) + H(\beta)(\beta^* - \beta) \]
- Estimasi parameter:
\[ \hat{\beta} \approx \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \]
9.3 Diagnostik Model GLM
Digunakan untuk mengevaluasi apakah model sudah sesuai dengan data.
- Uji formal
- Plot antara nilai prediksi vs aktual
Statistik Devians
- Mengukur apakah ada model lain yang lebih baik.
- Devians besar → model kurang cocok.
- Rumus:
\[ D = 2 \sum \left[ y_i \log\left( \frac{y_i}{\hat{\mu}_i} \right) - (y_i - \hat{\mu}_i) \right] \]
- Devians membandingkan model dengan saturated model.
- Devians kecil → model lebih cocok.
Statistik Chi-Kuadrat Pearson
- Mengukur apakah model lebih baik daripada tidak memakai model sama sekali:
\[ X^2 = \sum \frac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i} \]
- Jika signifikan → model memberikan informasi lebih baik.
Catatan
- Untuk data yang dikelompokkan, statistik mengikuti distribusi Chi-Square.
- Untuk data tidak dikelompokkan, tidak berlaku.
- Devians minimum dari MLE → alat evaluasi model.
Analisis Residual
- Residual adalah selisih antara nilai observasi dan prediksi.
- Berguna untuk menguji penyimpangan sistematis.
- Dapat digunakan untuk menguji asumsi model.
9.4 Estimasi dan Inferensi Regresi Logistik
Regresi logistik memodelkan probabilitas dari variabel respons biner (0/1) terhadap satu atau lebih variabel prediktor, dengan menggunakan MLE.
Fungsi Model Logistik
\[ \pi(x) = \frac{\exp(\beta_0 + \beta_1 x)}{1 + \exp(\beta_0 + \beta_1 x)} \]
Fungsi Log-Likelihood (n observasi):
\[ \ell(\beta) = \sum_{i=1}^n [y_i \log(\pi_i) + (1 - y_i)\log(1 - \pi_i)] \]
Estimasi dengan Newton-Raphson
Model regresi logistik:
\[ \pi_i = \frac{1}{1 + \exp(-x_i^T \beta)} \]
Log-likelihood:
\[ \ell(\beta) = \sum_{i=1}^n [y_i \log(\pi_i) + (1 - y_i)\log(1 - \pi_i)] \] di mana: - \(y_i\) adalah variabel respon biner (0 atau 1), - \(\pi_i = P(y_i = 1 | x_i) = \frac{1}{1 + e^{-x_i^\top \beta}}\) adalah probabilitas dari keberhasilan berdasarkan model logistik.
Langkah-Langkah Metode Newton-Raphson
Metode Newton-Raphson digunakan untuk mengestimasi parameter \(\beta\) secara iteratif dengan memperbarui nilai-nilai berdasarkan informasi dari turunan pertama dan kedua.
- Turunan Pertama (Score Function)
Turunan pertama dari log-likelihood terhadap parameter \(\beta\) adalah:
\[ U(\beta) = \frac{\partial \ell(\beta)}{\partial \beta} = X^\top (y - \pi) \]
di mana: - \(X\) adalah matriks desain (predictor), - \(y\) adalah vektor hasil observasi, - \(\pi\) adalah vektor probabilitas yang dihasilkan dari fungsi logit.
- Turunan Kedua (Hessian Matrix)
Hessian (turunan kedua) digunakan untuk mendapatkan arah pembaruan parameter dan didefinisikan sebagai:
\[ H(\beta) = -X^\top W X \]
dengan: - \(W\) adalah matriks diagonal dengan elemen \(\pi_i (1 - \pi_i)\).
- Iterasi Newton-Raphson
Pembaruan parameter dilakukan menggunakan formula berikut:
\[ \beta^{(t+1)} = \beta^{(t)} + (X^\top W^{(t)} X)^{-1} X^\top (y - \pi^{(t)}) \]
Proses ini diulang sampai konvergen, yaitu sampai perubahan nilai \(\beta\) cukup kecil.
Metode Newton-Raphson sangat penting dalam pendekatan numerik untuk mendapatkan estimasi parameter maksimum likelihood, terutama dalam regresi logistik yang tidak memiliki solusi eksplisit.
Perhitungan dengan R
# Data dalam format long (respon dikodekan)
data <- data.frame(
Race = rep(c("White", "White", "Black", "Black"), each = 3),
Gender = rep(c("Female", "Male"), each = 3, times = 2),
Belief = rep(c("Yes", "Undecided", "No"), times = 4),
Count = c(371, 49, 74,
250, 45, 71,
64, 9, 15,
25, 5, 13)
)
# Variabel dummy
data$Race_black <- ifelse(data$Race == "Black", 1, 0)
data$Gender_male <- ifelse(data$Gender == "Male", 1, 0)
data$Belief <- factor(data$Belief, levels = c("Yes", "Undecided", "No"))
library(nnet)
# Multinomial regression (Yes sebagai baseline)
model_mnlogit <- multinom(Belief ~ Race_black + Gender_male, weights = Count, data = data)## # weights: 12 (6 variable)
## initial value 1088.724778
## iter 10 value 773.726582
## final value 773.726510
## convergedsummary(model_mnlogit)## Call:
## multinom(formula = Belief ~ Race_black + Gender_male, data = data,
## weights = Count)
##
## Coefficients:
## (Intercept) Race_black Gender_male
## Undecided -2.025347 0.07079525 0.3134942
## No -1.643377 0.34177884 0.4185514
##
## Std. Errors:
## (Intercept) Race_black Gender_male
## Undecided 0.1473695 0.3091933 0.2083283
## No 0.1240866 0.2370373 0.1712550
##
## Residual Deviance: 1547.453
## AIC: 1559.453library(boot)
# Fungsi untuk bootstrap MLE
boot_fn <- function(data, indices) {
d <- data[indices, ]
fit <- multinom(Belief ~ Race_black + Gender_male, weights = Count, data = d, trace = FALSE)
coef(fit)
}
set.seed(123)
results <- boot(data = data, statistic = boot_fn, R = 100)## Warning in multinom(Belief ~ Race_black + Gender_male, weights = Count, : group
## 'Undecided' is empty## Warning in multinom(Belief ~ Race_black + Gender_male, weights = Count, : group
## 'Yes' is emptyresults##
## ORDINARY NONPARAMETRIC BOOTSTRAP
##
##
## Call:
## boot(data = data, statistic = boot_fn, R = 100)
##
##
## Bootstrap Statistics :
## original bias std. error
## t1* -2.02534716 1.1344441 9.156306
## t2* -1.64337735 0.7607659 10.251144
## t3* 0.07079525 -2.1594080 11.003624
## t4* 0.34177884 -2.0576244 10.358025
## t5* 0.31349418 0.2630050 12.017501
## t6* 0.41855145 0.4766176 11.707053# Uji Wald
z_values <- summary(model_mnlogit)$coefficients / summary(model_mnlogit)$standard.errors
p_values <- 2 * (1 - pnorm(abs(z_values)))
p_values## (Intercept) Race_black Gender_male
## Undecided 0 0.8188941 0.13237332
## No 0 0.1493368 0.01452445# Model tanpa prediktor (intercept-only)
model_null <- multinom(Belief ~ 1, weights = Count, data = data)## # weights: 6 (2 variable)
## initial value 1088.724778
## final value 778.098356
## converged# Likelihood Ratio Test
LR_stat <- 2 * (logLik(model_mnlogit) - logLik(model_null))
p_value_LR <- 1 - pchisq(LR_stat, df = length(coef(model_mnlogit)))
LR_stat## 'log Lik.' 8.743693 (df=6)p_value_LR## 'log Lik.' 0.188514 (df=6)# Akaike Information Criterion
AIC(model_mnlogit)## [1] 1559.453# Pseudo R-squared
ll_full <- as.numeric(logLik(model_mnlogit))
ll_null <- as.numeric(logLik(model_null))
pseudo_r2 <- 1 - (ll_full / ll_null)
pseudo_r2## [1] 0.005618629Interpretasi:
Berdasarkan hasil regresi logistik multinomial, dapat disimpulkan bahwa jenis kelamin memiliki pengaruh signifikan terhadap pilihan keyakinan terhadap kehidupan setelah mati, khususnya pada respon “No”. Individu berjenis kelamin laki-laki cenderung memiliki kemungkinan yang lebih tinggi untuk menyatakan tidak percaya pada kehidupan setelah mati dibandingkan perempuan, dengan nilai p yang signifikan (p = 0.014). Sebaliknya, variabel ras (Black vs White) tidak menunjukkan pengaruh yang signifikan terhadap pilihan keyakinan, baik pada respon “Undecided” maupun “No” dibandingkan dengan “Yes” sebagai kategori dasar.
Meskipun model berhasil diestimasi dan dikonfirmasi secara statistik melalui uji Wald dan Likelihood Ratio Test, kebaikan model secara keseluruhan tergolong rendah. Nilai pseudo R² McFadden yang sangat kecil (0.0056) menunjukkan bahwa model hanya menjelaskan kurang dari 1% variasi dalam data. Selain itu, hasil uji Likelihood Ratio tidak menunjukkan adanya peningkatan model yang signifikan dibandingkan model tanpa prediktor (intercept-only). Hal ini mengindikasikan bahwa model dengan prediktor ras dan jenis kelamin belum mampu menjelaskan preferensi kepercayaan terhadap kehidupan setelah mati secara memadai, kemungkinan disebabkan oleh ketidakseimbangan jumlah observasi antar kelompok atau kurangnya variabel penting lain dalam model.
9.5 Detail Metode Estimasi dan Inferensi Regresi Poisson
Regresi Poisson digunakan untuk menganalisis data bertipe hitung (jumlah kejadian) di mana variabel respon mengikuti distribusi Poisson. Estimasi parameter dilakukan menggunakan pendekatan Maximum Likelihood Estimation (MLE) dan inferensi statistik dilakukan melalui uji Wald dan Likelihood Ratio Test.
Model Regresi Poisson
Distribusi Poisson untuk sebuah kejadian \(Y_i\) didefinisikan sebagai:
\[ P(Y_i = y_i) = \frac{e^{-\lambda_i} \lambda_i^{y_i}}{y_i!} \]
Model regresi Poisson mengasumsikan bahwa:
\[ \log(\lambda_i) = x_i^\top \beta \]
Sehingga, nilai ekspektasi dari \(Y_i\) adalah:
\[ \lambda_i = \exp(x_i^\top \beta) \] Estimasi Parameter dengan MLE
Fungsi log-likelihood dari model regresi Poisson dapat dituliskan sebagai:
\[ \ell(\beta) = \sum_{i=1}^n \left[ y_i \log(\lambda_i) - \lambda_i - \log(y_i!) \right] \]
Dengan \(\lambda_i = \exp(x_i^\top \beta)\)
Estimasi Menggunakan Iterasi IRLS
Estimasi parameter \(\beta\) dilakukan secara iteratif menggunakan algoritma Iteratively Reweighted Least Squares (IRLS).
Tahap 1: Definisi Model
\[ \log(\lambda_i) = x_i^\top \beta \quad \Rightarrow \quad \lambda_i = \exp(x_i^\top \beta) \]
Tahap 2: Maksimisasi Log-Likelihood
\[ \ell(\beta) = \sum_{i=1}^n \left[ y_i \log(\lambda_i) - \lambda_i - \log(y_i!) \right] \] Tahap 3: Formulasi Iteratif IRLS
\[ \beta^{(t+1)} = \left( X^\top W^{(t)} X \right)^{-1} X^\top W^{(t)} z^{(t)} \]
dengan:
- \(W = \text{diag}(\lambda_i)\)
- \(z = \eta + \frac{y - \lambda}{\lambda}\), yaitu variabel kerja
- \(\eta = \log(\lambda) = X^\top \beta\)
Metode IRLS sangat berguna dalam memperoleh estimasi parameter regresi Poisson secara numerik ketika tidak tersedia solusi analitik.
Perhitungan dengan R
# Subset data hanya "Yes" dan "No"
data_bin <- subset(data, Belief != "Undecided")
data_bin$Belief_bin <- ifelse(data_bin$Belief == "Yes", 1, 0) # 1 = Yes, 0 = No
#Model GLM
model_bin <- glm(Belief_bin ~ Race_black + Gender_male,
data = data_bin,
weights = Count,
family = binomial)
summary(model_bin)##
## Call:
## glm(formula = Belief_bin ~ Race_black + Gender_male, family = binomial,
## data = data_bin, weights = Count)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.6431 0.1241 13.244 <2e-16 ***
## Race_black -0.3409 0.2371 -1.438 0.1505
## Gender_male -0.4182 0.1713 -2.442 0.0146 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 873.64 on 7 degrees of freedom
## Residual deviance: 866.16 on 5 degrees of freedom
## AIC: 872.16
##
## Number of Fisher Scoring iterations: 5# Uji Wald
coefs <- coef(summary(model_bin))
wald_z <- coefs[, "Estimate"] / coefs[, "Std. Error"]
wald_p <- 2 * (1 - pnorm(abs(wald_z)))
wald_p## (Intercept) Race_black Gender_male
## 0.00000000 0.15046944 0.01461569# Uji Likehood Ratio
model_null_bin <- glm(Belief_bin ~ 1, data = data_bin, weights = Count, family = binomial)
lr_stat_bin <- 2 * (logLik(model_bin) - logLik(model_null_bin))
p_val_lr_bin <- 1 - pchisq(lr_stat_bin, df = length(coef(model_bin)) - 1)
lr_stat_bin## 'log Lik.' 7.481653 (df=3)p_val_lr_bin## 'log Lik.' 0.02373447 (df=3)#Evaluasi Model
AIC(model_bin)## [1] 872.1555ll_full_bin <- as.numeric(logLik(model_bin))
ll_null_bin <- as.numeric(logLik(model_null_bin))
pseudo_r2_bin <- 1 - (ll_full_bin / ll_null_bin)
pseudo_r2_bin## [1] 0.0085638Interpretasi:
Berdasarkan hasil regresi logistik biner untuk membedakan antara respon “Yes” dan “No” terhadap keyakinan akan kehidupan setelah mati, model berhasil dikonvergensi dalam 5 iterasi menggunakan metode IRLS. Koefisien intercept sebesar 1.6431 signifikan pada tingkat signifikansi 0,001, menunjukkan bahwa kelompok referensi (White Female) memiliki peluang yang lebih tinggi untuk menjawab “Yes”. Variabel Race_black tidak signifikan (p = 0.1505), menunjukkan bahwa tidak terdapat perbedaan yang berarti antara ras kulit putih dan kulit hitam dalam hal keyakinan tersebut setelah mengendalikan gender. Sebaliknya, Gender_male memiliki pengaruh negatif dan signifikan (p = 0.0146), yang mengindikasikan bahwa laki-laki cenderung lebih kecil kemungkinannya untuk menjawab “Yes” dibandingkan perempuan, dengan asumsi ras yang sama.
Hasil uji Likelihood Ratio menunjukkan bahwa model dengan prediktor (ras dan gender) secara keseluruhan signifikan lebih baik dibanding model null (hanya intercept), dengan statistik LR sebesar 7.48 dan nilai p = 0.0237. Meskipun demikian, nilai pseudo R² McFadden relatif kecil (0.0086), yang menandakan bahwa variabel yang digunakan dalam model hanya menjelaskan kurang dari 1% variasi dalam data. Hal ini menunjukkan bahwa meskipun ada pengaruh yang signifikan dari gender, secara keseluruhan model ini belum cukup kuat untuk menjelaskan keputusan individu dalam mempercayai kehidupan setelah mati, dan mungkin diperlukan tambahan variabel lain untuk meningkatkan daya jelaskan model.
Bab X Regresi Logistik dengan Prediktor Nominal, Ordinal, dan Rasio
Regresi logistik digunakan untuk memodelkan hubungan antara variabel respons biner dengan satu atau lebih prediktor. Tiga jenis skala pengukuran prediktor dapat digunakan:
- Nominal: Kategori tanpa urutan (misalnya: jenis kelamin, ras).
- Ordinal: Kategori dengan urutan tetapi tanpa jarak numerik yang pasti (misalnya: tingkat pendidikan).
- Rasio: Data kuantitatif dengan nol absolut dan rasio bermakna (misalnya: pendapatan).
Regresi logistik biner menghitung log-odds (logaritma dari peluang keberhasilan) sebagai fungsi linear dari prediktor.
10.1 Simulasi Data
Simulasi data dengan 500 observasi menggunakan tiga variabel prediktor: sektor pekerjaan (nominal), tingkat kepuasan (ordinal), dan lama bekerja (rasio).
library(tibble)
library(knitr)
# Simulasi variabel
n <- 500
# Variabel prediktor
set.seed(2025)
sektor <- sample(c("Pemerintah", "Swasta", "Wirausaha"), n, replace = TRUE)
kepuasan <- sample(c("Rendah", "Sedang", "Tinggi"), n, replace = TRUE, prob = c(0.3, 0.5, 0.2))
lama_kerja <- rnorm(n, mean = 7, sd = 3)
# Model logit untuk simulasi loyalitas
logit_p <- -1 + 0.4 * (sektor == "Swasta") +
0.8 * as.numeric(factor(kepuasan, levels = c("Rendah", "Sedang", "Tinggi"), ordered = TRUE)) +
0.1 * lama_kerja
p <- 1 / (1 + exp(-logit_p))
# Variabel respons
loyal <- rbinom(n, 1, p)
# Dataset
sim_data <- tibble(loyal, sektor, kepuasan, lama_kerja)
head(sim_data)Interpretasi: Dataset ini berisi 500 observasi
dengan variabel loyal, sektor,
kepuasan, dan lama_kerja. Variabel
loyal disimulasikan menggunakan fungsi logistik yang
dipengaruhi oleh ketiga prediktor tersebut.
10.2 Eksplorasi Data
library(dplyr)
sim_data %>%
group_by(loyal) %>%
summarise(
Jumlah = n(),
Rata2_LamaKerja = mean(lama_kerja)
)Interpretasi: Distribusi jumlah loyal dan tidak loyal, serta rata-rata lama kerja pada masing-masing grup.
10.3 Perlakuan Variabel Ordinal
10.3.1 Treat sebagai Nominal (Dummy)
sim_data_nominal <- sim_data %>%
mutate(
kepuasan = factor(kepuasan, levels = c("Rendah", "Sedang", "Tinggi"))
)
model_nominal <- glm(loyal ~ sektor + kepuasan + lama_kerja,
data = sim_data_nominal,
family = binomial)
summary(model_nominal)##
## Call:
## glm(formula = loyal ~ sektor + kepuasan + lama_kerja, family = binomial,
## data = sim_data_nominal)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.24786 0.35445 -0.699 0.4844
## sektorSwasta 0.42545 0.28843 1.475 0.1402
## sektorWirausaha 0.29545 0.27396 1.078 0.2808
## kepuasanSedang 0.96515 0.23643 4.082 4.46e-05 ***
## kepuasanTinggi 2.53828 0.54003 4.700 2.60e-06 ***
## lama_kerja 0.08164 0.03855 2.118 0.0342 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 516.59 on 499 degrees of freedom
## Residual deviance: 467.96 on 494 degrees of freedom
## AIC: 479.96
##
## Number of Fisher Scoring iterations: 5Interpretasi Koefisien:
- Intercept:
Koefisien = -0.24786 → baseline log-odds loyal saat sektor = Pemerintah, kepuasan = Rendah, lama kerja = 0. Tidak signifikan (p = 0.4844).
- Sektor:
sektorSwasta → Koefisien = +0.42545 dan odds ratio ≈ 1.53. Tidak signifikan (p = 0.1402). Artinya peluang loyal sedikit lebih tinggi di Swasta dibanding Pemerintah, tetapi tidak signifikan.
sektorWirausaha → Koefisien = +0.29545 dan odds ratio ≈ 1.34. Tidak signifikan (p = 0.2808). Artinya peluang loyal sedikit lebih tinggi di Wirausaha dibanding Pemerintah, tetapi tidak signifikan.
- Kepuasan:
kepuasanSedang → Koefisien = +0.96515 dan odds ratio ≈ 2.63. Sangat signifikan (p < 0.001). Peluang loyal 2.63 kali lebih besar dibanding kepuasan Rendah.
kepuasanTinggi → Koefisien = +2.53828 → odds ratio ≈ 12.64. Sangat signifikan (p < 0.001). Peluang loyal ~12.64 kali lebih besar dibanding kepuasan Rendah.
- Lama Kerja:
Koefisien = +0.08164 dan odds ratio ≈ 1.085. Signifikan (p = 0.0342). Setiap tambahan 1 tahun kerja, peluang loyal meningkat sekitar 8.5%.
Kesimpulan
Faktor signifikan yang mempengaruhi loyalitas, yaitu Kepuasan (efek sangat besar, terutama untuk kepuasan tinggi) dan Lama kerja (efek moderat). Sedangkan, faktor yang tidak signifikan, yaitu Sektor pekerjaan (Swasta, Wirausaha vs Pemerintah)
10.3.2 Treat sebagai Rasio (Numeric Rank)
sim_data_score <- sim_data %>%
mutate(
skor_kepuasan = recode(kepuasan,
"Rendah" = 1,
"Sedang" = 2,
"Tinggi" = 3)
)
model_score <- glm(loyal ~ sektor + skor_kepuasan + lama_kerja,
data = sim_data_score,
family = binomial)
summary(model_score)##
## Call:
## glm(formula = loyal ~ sektor + skor_kepuasan + lama_kerja, family = binomial,
## data = sim_data_score)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.40327 0.45355 -3.094 0.00197 **
## sektorSwasta 0.41364 0.28849 1.434 0.15163
## sektorWirausaha 0.29225 0.27446 1.065 0.28696
## skor_kepuasan 1.10697 0.18719 5.914 3.35e-09 ***
## lama_kerja 0.08232 0.03850 2.138 0.03250 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 516.59 on 499 degrees of freedom
## Residual deviance: 468.96 on 495 degrees of freedom
## AIC: 478.96
##
## Number of Fisher Scoring iterations: 4Interpretasi Koefisien:
- Intercept:
Koefisien = undefined → baseline log-odds loyal saat sektor = Pemerintah, kepuasan = 0, lama kerja = 0. Pada kondisi baseline, peluang loyal sangat rendah. P-value signifikan (p = 0.00197), artinya intercept berbeda signifikan dari nol.
- Sektor:
sektorSwasta → Koefisien = +0.41364 dan odds ratio ≈ exp(0.41364) ≈ 1.512. Tidak signifikan (p = 0.15163). Artinya peluang loyal di sektor Swasta 1.51 kali lebih besar dibanding Pemerintah, namun perbedaan ini tidak signifikan.
sektorWirausaha → Koefisien = +0.29225 dan odds ratio ≈ exp(0.29225) ≈ 1.339. Tidak signifikan (p = 0.28696). Artinya Peluang loyal di sektor Wirausaha 1.34 kali lebih besar dibanding Pemerintah, tetapi perbedaan ini tidak signifikan.
- Kepuasan:
Koefisien = +1.10697 → odds ratio ≈ exp(1.10697) ≈ 3.025. Sangat signifikan (p < 0.001). Artinya setiap peningkatan 1 poin skor_kepuasan meningkatkan peluang loyal ~3 kali lipat. Ini menunjukkan skor_kepuasan adalah prediktor kuat loyalitas.
- Lama Kerja:
Koefisien = +0.08232 dan odds ratio ≈ exp(0.08232) ≈ 1.086. Signifikan (p = 0.03250). Artinya setiap peningkatan 1 tahun lama kerja meningkatkan peluang loyal sekitar 8.6%.
Kesimpulan
Faktor signifikan yang mempengaruhi loyalitas, yaitu Skor kepuasan (semakin tinggi skor, semakin besar peluang loyal) dan Lama kerja (semakin lama kerja, semakin besar peluang loyal). Sedangkan, faktor yang tidak signifikan, yaitu Sektor pekerjaan (Swasta, Wirausaha vs Pemerintah tidak berpengaruh signifikan terhadap loyalitas).
10.3.3 Perbandingan Model
{AIC(model_nominal, model_score)}Interpretasi:
Nilai AIC digunakan untuk memilih model terbaik, di mana nilai lebih
kecil menunjukkan model yang lebih baik. Jika model_score
memiliki nilai AIC lebih kecil dari model_nominal, maka
perlakuan ordinal sebagai numerik (skor) menghasilkan model yang lebih
efisien.
10.3.4 Goodness of Fit (GOF) Model
# Model null (tanpa prediktor)
nullmod <- glm(loyal ~ 1, data = sim_data, family = binomial)
# McFadden's pseudo-R²
r2_nominal <- 1 - (logLik(model_nominal) / logLik(nullmod))
r2_score <- 1 - (logLik(model_score) / logLik(nullmod))
# Hasil
list(
McFadden_R2_Nominal = r2_nominal,
McFadden_R2_Score = r2_score
)## $McFadden_R2_Nominal
## 'log Lik.' 0.0941363 (df=6)
##
## $McFadden_R2_Score
## 'log Lik.' 0.09221591 (df=5)Interpretasi:
Kedua model (Nominal & Score) punya nilai McFadden R² < 0.1 → menunjukkan bahwa model belum optimal dalam menjelaskan variabel loyalitas. Namun, meskipun R² kecil, variabel skor kepuasan dan lama kerja tetap signifikan, sehingga ada kontribusi bermakna.
10.4 Visualisasi Prediksi
library(ggplot2)
# Visualisasi untuk model numerik
sim_data_score$predicted_prob <- predict(model_score, type = "response")
ggplot(sim_data_score, aes(x = lama_kerja, y = predicted_prob, color = sektor)) +
geom_point(alpha = 0.5) +
facet_wrap(~kepuasan) +
labs(title = "Prediksi Probabilitas Loyalitas (Model Numerik)",
y = "Probabilitas", x = "Lama Bekerja (tahun)") +
theme_minimal()
# Visualisasi untuk model nominal
sim_data_nominal$predicted_prob <- predict(model_nominal, type = "response")
ggplot(sim_data_nominal, aes(x = lama_kerja, y = predicted_prob, color = sektor)) +
geom_point(alpha = 0.5) +
facet_wrap(~kepuasan) +
labs(title = "Prediksi Probabilitas Loyalitas (Model Nominal)",
y = "Probabilitas", x = "Lama Bekerja (tahun)") +
theme_minimal()
Interpretasi:
Plot menunjukkan hubungan antara lama bekerja dan probabilitas loyalitas, dikelompokkan berdasarkan sektor dan kepuasan. Perbedaan pola antara model nominal dan numerik dapat membantu menentukan perlakuan variabel yang paling sesuai.
Bab XI Pemilihan Model Regresi Logistik dan Evaluasi
Dalam analisis regresi logistik, pemilihan model merupakan langkah penting untuk memperoleh model yang mampu memprediksi probabilitas terjadinya suatu peristiwa (respon biner) dengan baik. Terdapat dua pendekatan utama yang dapat digunakan, yaitu pendekatan Confirmatory dan Exploratory.
11.1 Membangun Model Regresi Logistik: Pendekatan Confirmatory dan Exploratory
11.1.1 Confirmatory (Pendekatan Konfirmatori)
Pendekatan ini diterapkan ketika peneliti telah memiliki landasan teori atau hipotesis yang jelas mengenai pengaruh atau hubungan antar variabel prediktor dan respon.
Ciri-ciri:
- Model dikembangkan berdasarkan teori yang ada atau temuan penelitian sebelumnya.
- Tujuan utama adalah menguji apakah pengaruh variabel sesuai dengan teori dan signifikan, bukan sekadar mencari model dengan performa terbaik.
- Umumnya, peneliti membangun model penuh terlebih dahulu, kemudian mengevaluasi apakah penambahan atau pengurangan variabel (misalnya interaksi) memberikan peningkatan model yang bermakna.
- Pengujian signifikansi dilakukan dengan membandingkan model dengan dan tanpa variabel tertentu menggunakan uji seperti Likelihood Ratio Test.
Contoh penggunaan:
Misalnya, sebuah teori menyatakan bahwa faktor x1
dan x2 mempengaruhi probabilitas seseorang membeli
produk.
Model regresi logistik dibangun dengan memasukkan variabel
x1 dan x2, kemudian diuji apakah
keduanya memberikan kontribusi signifikan terhadap prediksi probabilitas
tersebut.
11.1.2 Exploratory (Pendekatan Eksploratori)
Pendekatan ini digunakan ketika peneliti belum memiliki teori yang pasti, atau bertujuan untuk mengeksplorasi hubungan potensial antar variabel.
Ciri-ciri:
- Model dikembangkan secara bertahap dengan tujuan menemukan kombinasi prediktor yang paling optimal.
- Pemilihan variabel dilakukan berdasarkan kriteria statistik, seperti AIC, deviance, atau log-likelihood.
Proses seleksi variabel dapat dilakukan melalui:
- Forward Selection: Memulai dari model kosong, menambahkan variabel satu per satu jika terbukti signifikan.
- Backward Elimination: Memulai dari model penuh, secara bertahap mengeluarkan variabel yang tidak signifikan.
- Stepwise Selection: Mengombinasikan forward dan backward; variabel dapat masuk dan keluar model secara dinamis.
Tujuan:
Menemukan model yang parsimonious — model yang sederhana namun tetap memiliki performa prediksi yang baik.
11.2 Pemilihan Model
Confirmatory (Model Full)
model_full <- glm(loyal ~ sektor + kepuasan + lama_kerja, data = sim_data, family = binomial)
summary(model_full)##
## Call:
## glm(formula = loyal ~ sektor + kepuasan + lama_kerja, family = binomial,
## data = sim_data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.24786 0.35445 -0.699 0.4844
## sektorSwasta 0.42545 0.28843 1.475 0.1402
## sektorWirausaha 0.29545 0.27396 1.078 0.2808
## kepuasanSedang 0.96515 0.23643 4.082 4.46e-05 ***
## kepuasanTinggi 2.53828 0.54003 4.700 2.60e-06 ***
## lama_kerja 0.08164 0.03855 2.118 0.0342 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 516.59 on 499 degrees of freedom
## Residual deviance: 467.96 on 494 degrees of freedom
## AIC: 479.96
##
## Number of Fisher Scoring iterations: 5Exploratory (Stepwise AIC)
null_model <- glm(loyal ~ 1, data = sim_data, family = binomial)
step_forward <- step(null_model, direction = "forward", scope = formula(model_full), trace = FALSE)
step_backward <- step(model_full, direction = "backward", trace = FALSE)
step_both <- step(null_model, direction = "both", scope = formula(model_full), trace = FALSE)
AIC(model_full, step_forward, step_backward, step_both)Interpretasi:
- Model Full: AIC = 479.96 dan McFadden R² = 0.0941, artinya model menjelaskan sekitar 9.4% variasi loyalitas (masih tergolong kecil).
- Model Stepwise: AIC = 478.96 dan McFadden R² = 0.0922, sedikit lebih sederhana dan lebih optimal. Model stepwise cenderung memilih model yang lebih parsimonious.
11.3 Evaluasi Model: ROC dan AUC
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, varpred_prob <- predict(step_both, type = "response")
roc_obj <- roc(sim_data$loyal, pred_prob)## Setting levels: control = 0, case = 1## Setting direction: controls < casesplot(roc_obj, main = "ROC Curve - Stepwise Both Model", col = "darkblue")
auc(roc_obj)## Area under the curve: 0.6996Interpretasi:
AUC = 0.6996 artinya model memiliki kemampuan klasifikasi yang sedang (fair). 0.7–0.8 = fair artinya model ini cukup berguna, tetapi belum sangat kuat.
11.4 Pseudo R-Squared
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(step_both)[c("McFadden", "r2ML", "r2CU")] # McFadden, Cox-Snell (ML), Nagelkerke (CU)## fitting null model for pseudo-r2## McFadden r2ML r2CU
## 0.08961512 0.08843225 0.13728950Interpretasi:
Model Nagelkerke (CU) yang paling optimal dibandingkan yang lain.
11.5 Tabel Klasifikasi dan Evaluasi
library(caret)## Warning: package 'caret' was built under R version 4.3.3## Loading required package: lattice##
## Attaching package: 'lattice'## The following object is masked from 'package:boot':
##
## melanoma##
## Attaching package: 'caret'## The following object is masked from 'package:purrr':
##
## liftpred_class <- ifelse(pred_prob >= 0.5, 1, 0)
conf_matrix <- confusionMatrix(factor(pred_class), factor(sim_data$loyal), positive = "1"); conf_matrix## Confusion Matrix and Statistics
##
## Reference
## Prediction 0 1
## 0 1 1
## 1 105 393
##
## Accuracy : 0.788
## 95% CI : (0.7495, 0.823)
## No Information Rate : 0.788
## P-Value [Acc > NIR] : 0.526
##
## Kappa : 0.0108
##
## Mcnemar's Test P-Value : <2e-16
##
## Sensitivity : 0.997462
## Specificity : 0.009434
## Pos Pred Value : 0.789157
## Neg Pred Value : 0.500000
## Prevalence : 0.788000
## Detection Rate : 0.786000
## Detection Prevalence : 0.996000
## Balanced Accuracy : 0.503448
##
## 'Positive' Class : 1
## conf_matrix$byClass[c("Sensitivity", "Specificity")]## Sensitivity Specificity
## 0.997461929 0.009433962Interpretasi:
Model sangat baik dalam memprediksi Loyal (Sensitivity tinggi). Akurasi terlihat tinggi (78.8%) hanya karena mayoritas data memang Loyal. Balanced Accuracy hanya 50.3% artinya secara adil untuk kedua kelas, performanya hampir sama dengan acak.
11.6 Metode Perbandingan Model dalam Regresi Logistik
# Model 1: loyal ~ sektor (satu prediktor)
model1 <- glm(loyal ~ sektor, data = sim_data, family = binomial)
# Model 2: loyal ~ sektor + kepuasan (dua prediktor)
model2 <- glm(loyal ~ sektor + kepuasan, data = sim_data, family = binomial)
# Model 3: loyal ~ sektor + kepuasan + lama_kerja (tiga prediktor)
model3 <- glm(loyal ~ sektor + kepuasan + lama_kerja, data = sim_data, family = binomial)
# Perbandingan AIC dan Deviance
model_comp <- data.frame(
Model = c("Model 1", "Model 2", "Model 3"),
AIC = c(AIC(model1), AIC(model2), AIC(model3)),
Deviance = c(deviance(model1), deviance(model2), deviance(model3))
) ; model_comp11.7 Likelihood-Ratio Test
anova(model1, model2, test = "LRT")anova(model2, model3, test = "LRT")11.8 Prinsip Parsimony
Dalam pemilihan model, model yang terlalu kompleks memang cenderung memiliki nilai AIC dan deviance yang lebih rendah. Namun,
Rumus AIC
\[ \text{AIC} = -2(\log L - k) = -2 \log L + 2k \]
AIC merupakan ukuran yang mengombinasikan kualitas fit model (melalui log-likelihood) dan tingkat kompleksitas model (melalui jumlah parameter \(k\)). Semakin kecil nilai AIC, semakin baik model tersebut secara keseluruhan. Namun, AIC memberikan penalti bagi model yang terlalu kompleks agar tidak sekadar mengejar likelihood tinggi.
Rumus Deviance
\[ D = -2 [\log L(\text{model}) - \log L(\text{model saturasi})] \]
Deviance mengukur seberapa jauh prediksi model saat ini dibandingkan dengan model sempurna (model saturasi). Semakin kecil nilai deviance, semakin baik model dalam merepresentasikan data aktual.
Rumus Likelihood-Ratio
\[ G^2 = -2 [\log L_0 - \log L_1] \]
Statistik Likelihood Ratio digunakan untuk menguji apakah penambahan variabel dalam model secara signifikan meningkatkan kecocokan model. Jika nilai \(G^2\) besar dan p-value kecil, maka model yang lebih kompleks memberikan fit yang lebih baik dibanding model yang lebih sederhana.
11.9 Evaluasi Tabel Klasifikasi dan Akurasi Model
library(caret)
pred_prob <- predict(model3, type = "response")
pred_class <- factor(ifelse(pred_prob >= 0.5, 1, 0))
sim_data$loyal <- factor(sim_data$loyal)
conf_matrix <- confusionMatrix(pred_class, sim_data$loyal, positive = "1")
conf_matrix## Confusion Matrix and Statistics
##
## Reference
## Prediction 0 1
## 0 2 0
## 1 104 394
##
## Accuracy : 0.792
## 95% CI : (0.7537, 0.8268)
## No Information Rate : 0.788
## P-Value [Acc > NIR] : 0.4389
##
## Kappa : 0.0294
##
## Mcnemar's Test P-Value : <2e-16
##
## Sensitivity : 1.00000
## Specificity : 0.01887
## Pos Pred Value : 0.79116
## Neg Pred Value : 1.00000
## Prevalence : 0.78800
## Detection Rate : 0.78800
## Detection Prevalence : 0.99600
## Balanced Accuracy : 0.50943
##
## 'Positive' Class : 1
## 11.9.1 Sensitivitas dan Spesifisitas
Sensitivitas Kemampuan model untuk secara benar mengidentifikasi kelas positif (True Positive Rate).
\[ \text{Sensitivity} = \frac{\text{TP}}{\text{TP} + \text{FN}} \]
Spesifisitas Kemampuan model untuk secara benar mengidentifikasi kelas negatif (True Negative Rate).
\[ \text{Specificity} = \frac{\text{TN}}{\text{TN} + \text{FP}} \]
conf_matrix$byClass[c("Sensitivity", "Specificity")]## Sensitivity Specificity
## 1.00000000 0.01886792Interpretasi:
- Akurasi (Accuracy): 79.2%. Artinya, model berhasil memprediksi benar sekitar 79.2% dari seluruh data. Namun, karena data sangat imbalanced (kelas 1 jauh lebih dominan), akurasi tinggi ini bisa menyesatkan.
- Sensitivitas (Sensitivity) = 100%. Nilai 1.00 berarti model tidak pernah gagal mendeteksi kelas positif — seluruh observasi yang benar-benar loyal diprediksi sebagai loyal (TP = 394, FN = 0).
- Spesifisitas (Specificity) ≈ 1.89%. Nilai sangat rendah (1.89%) artinya model hampir selalu salah memprediksi kelas 0 — dari 106 observasi loyal = 0, hanya 2 yang berhasil diprediksi benar (TN = 2), sisanya (104) diprediksi sebagai loyal = 1 (False Positive).
- Model sangat bagus dalam mendeteksi pelanggan loyal (positif). Sedangkan, Model sangat buruk dalam mendeteksi pelanggan tidak loyal (negatif) → Spesifisitas hampir nol. Artinya, jika digunakan untuk memprediksi siapa yang tidak loyal, model ini tidak bisa diandalkan.
11.10 ROC (Receiver Operating Characteristic)
Kurva ROC merupakan alat visual yang digunakan untuk mengevaluasi kinerja model klasifikasi biner. Kurva ini menggambarkan perbandingan antara True Positive Rate (Sensitivity) dengan False Positive Rate (1 - Specificity) pada berbagai ambang klasifikasi.
11.10.1 Definisi
- Sumbu Y (Sensitivity / True Positive Rate)
Kemampuan model untuk mengidentifikasi kelas positif dengan benar.
\[ \text{Sensitivity} = \frac{TP}{TP + FN} \]
- Sumbu X (1 - Specificity / False Positive
Rate)
Proporsi kasus negatif yang secara keliru diklasifikasikan sebagai positif.
\[ 1 - \text{Specificity} = \frac{FP}{FP + TN} \]
Garis Diagonal
Menunjukkan performa acak (random guess).Kurva yang mendekati sudut kiri atas
Mengindikasikan performa klasifikasi yang lebih baik.
11.10.2 Cut-off dan Pergerakan Kurva
- Cut-off rendah → model lebih liberal:
- Sensitivitas meningkat
- Spesifisitas menurun
- Cut-off tinggi → model lebih konservatif:
- Sensitivitas menurun
- Spesifisitas meningkat
11.10.3 ROC Ideal
Kurva ROC yang ideal ditandai oleh:
- Garis vertikal tajam menuju sensitivitas = 1
- Garis horizontal mendekati spesifisitas = 1
- Area Under the Curve (AUC) mendekati 1
11.10.4 Interpretasi Area Under the Curve (AUC)
- AUC = 0.5 → model sama saja dengan tebakan acak.
- AUC > 0.7 → model cukup baik.
- AUC > 0.9 → model sangat baik.
AUC juga dikenal sebagai concordance index, yaitu probabilitas bahwa model memberikan skor prediksi lebih tinggi untuk kasus positif dibandingkan kasus negatif.
11.10.5 Kegunaan Kurva ROC
- Untuk membandingkan performa antar model klasifikasi.
- Untuk memilih threshold (cut-off) optimal sesuai dengan kebutuhan
aplikasi.
Contoh: apakah lebih penting menghindari false negative atau false positive?
11.10.6 Simulasi Pemilihan Threshold Optimal
thresholds <- seq(0.1, 0.9, by = 0.05)
results <- data.frame(Threshold = thresholds)
obs <- factor(sim_data$loyal, levels = c("0", "1"))
results$Sensitivity <- sapply(thresholds, function(t) {
pred_class <- factor(ifelse(pred_prob >= t, "1", "0"), levels = c("0", "1"))
cm <- table(Pred = pred_class, Obs = obs)
TP <- cm["1", "1"]
FN <- cm["0", "1"]
TP / (TP + FN)
})
results$Specificity <- sapply(thresholds, function(t) {
pred_class <- factor(ifelse(pred_prob >= t, "1", "0"), levels = c("0", "1"))
cm <- table(Pred = pred_class, Obs = obs)
TN <- cm["0", "0"]
FP <- cm["1", "0"]
TN / (TN + FP)
})
print(results)## Threshold Sensitivity Specificity
## 1 0.10 1.0000000 0.00000000
## 2 0.15 1.0000000 0.00000000
## 3 0.20 1.0000000 0.00000000
## 4 0.25 1.0000000 0.00000000
## 5 0.30 1.0000000 0.00000000
## 6 0.35 1.0000000 0.00000000
## 7 0.40 1.0000000 0.00000000
## 8 0.45 1.0000000 0.00000000
## 9 0.50 1.0000000 0.01886792
## 10 0.55 0.9670051 0.10377358
## 11 0.60 0.9213198 0.20754717
## 12 0.65 0.8756345 0.30188679
## 13 0.70 0.7994924 0.43396226
## 14 0.75 0.7335025 0.52830189
## 15 0.80 0.6065990 0.70754717
## 16 0.85 0.3578680 0.89622642
## 17 0.90 0.2309645 0.96226415Interpretasi:
Berdasarkan analisis Sensitivity dan Specificity pada berbagai threshold, cut-off optimal diperoleh pada threshold sekitar 0.80, di mana kombinasi Sensitivity + Specificity mencapai maksimum. Jika tujuan aplikasi adalah meminimalkan false negative, maka threshold sekitar 0.65–0.70 dapat dipertimbangkan agar sensitivitas tetap tinggi. Namun, karena data ini tampak imbalance (kelas 1 lebih dominan), visualisasi tambahan dengan Precision-Recall Curve disarankan.
11.11 Precision-Recall Curve (PR Curve)
Kurva Precision-Recall (PR) adalah alat visual untuk mengevaluasi kinerja model klasifikasi, terutama sangat berguna ketika menghadapi data yang tidak seimbang (class imbalance).
11.11.1 Definisi
- Precision (Presisi): Proporsi prediksi positif yang benar-benar positif.
\[ Precision = \frac{TP}{TP + FP} \]
- Recall (Sensitivitas): Proporsi kasus positif yang berhasil diprediksi sebagai positif.
\[ Recall = \frac{TP}{TP + FN} \]
11.11.2 Interpretasi
- Kurva PR menunjukkan bagaimana nilai presisi berubah seiring dengan meningkatnya recall.
- Idealnya, kita menginginkan precision dan recall yang tinggi secara bersamaan, walaupun seringkali ada trade-off antara keduanya.
- Model yang berkinerja baik akan menghasilkan kurva PR yang melengkung mendekati pojok kanan atas.
11.11.3 Area Under PR Curve (AUPRC)
- Luas kurva PR yang mendekati 1 menunjukkan model yang sangat baik.
- Baseline AUPRC setara dengan prevalensi kelas positif di dalam data.
11.11.4 PR Curve vs ROC Curve
| Aspek | ROC Curve | Precision-Recall Curve |
|---|---|---|
| Fokus | Semua kelas | Hanya kelas positif |
| Cocok digunakan untuk | Data seimbang | Data tidak seimbang |
| Sumbu Y | Sensitivitas (Recall) | Precision |
| Sumbu X | 1 - Spesifisitas |
11.11.5 Visualisasi PR Curve
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##
## Attaching package: 'rlang'## The following objects are masked from 'package:purrr':
##
## %@%, flatten, flatten_chr, flatten_dbl, flatten_int, flatten_lgl,
## flatten_raw, invoke, spliceprob <- pred_prob
label <- as.numeric(as.character(sim_data$loyal))
pr <- pr.curve(
scores.class0 = prob[label == 1],
scores.class1 = prob[label == 0],
curve = TRUE
)
plot(pr)
Interpretasi:
Kurva precision-recall yang tinggi pada recall rendah menunjukkan bahwa prediksi loyalitas yang paling yakin dari model sangat akurat, sementara trade-off precision yang menurun saat recall meningkat mencerminkan upaya model untuk menangkap lebih banyak individu loyal dengan risiko meningkatnya kesalahan.
11.12 Pseudo R-squared pada Regresi Logistik
Pada regresi logistik, nilai R-squared konvensional tidak dapat digunakan seperti dalam regresi linear. Oleh karena itu, digunakan ukuran alternatif yang dikenal sebagai pseudo R-squared. Dua versi umum dari pseudo R-squared adalah:
- \(R^2_{\text{Cox and Snell}}\)
- \(R^2_{\text{McFadden}}\)
model <- model3
model_null <- glm(loyal ~ 1, data = sim_data, family = binomial)Likelihood
Berikut adalah rumus untuk masing-masing pseudo R-squared:
\[ R^2_{\text{Cox and Snell}} = 1 - \left(\frac{L_0}{L_M}\right)^{2/n} \]
\[ R^2_{\text{McFadden}} = 1 - \frac{\log L_M}{\log L_0} \]
Dengan:
- \(L_0\): nilai likelihood dari model nol (tanpa prediktor)
- \(L_M\): nilai likelihood dari model penuh (dengan prediktor)
- \(n\): jumlah observasi
11.12.1 Perhitungan Manual
logL0 <- logLik(model_null)
logLM <- logLik(model)
L0 <- exp(logL0)
LM <- exp(logLM)
n <- nobs(model)
cox_snell <- 1 - (L0 / LM)^(2 / n)
mcfadden <- 1 - (as.numeric(logLM) / as.numeric(logL0))
r2 <- data.frame(
R2_Cox_Snell = cox_snell,
R2_McFadden = mcfadden
)
r211.12.2 Perhitungan Otomatis dengan Package Tambahan
# Pseudo R2 dengan pscl
if (!require(pscl)) install.packages("pscl", dependencies = TRUE)
library(pscl)
pR2(model)## fitting null model for pseudo-r2## llh llhNull G2 McFadden r2ML
## -233.98209894 -258.29724697 48.63029607 0.09413630 0.09268047
## r2CU
## 0.14388478# Nagelkerke R² dengan rcompanion
if (!require(DescTools)) install.packages("DescTools", dependencies = TRUE)## Loading required package: DescTools## Warning: package 'DescTools' was built under R version 4.3.3##
## Attaching package: 'DescTools'## The following objects are masked from 'package:caret':
##
## MAE, RMSElibrary(DescTools)
PseudoR2(model3, which = "Nagelkerke")## Nagelkerke
## 0.1438848# Semua jenis Pseudo R² dengan DescTools
if (!require(DescTools)) install.packages("DescTools", dependencies = TRUE)
library(DescTools)
PseudoR2(model, which = "all")## McFadden McFaddenAdj CoxSnell Nagelkerke AldrichNelson
## 0.09413630 0.07090725 0.09268047 0.14388478 0.08863947
## VeallZimmermann Efron McKelveyZavoina Tjur AIC
## 0.17443158 0.08541782 0.20301960 0.08906180 479.96419787
## BIC logLik logLik0 G2
## 505.25184646 -233.98209894 -258.29724697 48.63029607Interpretasi:
- Nagelkerke R² = 0.14 → sekitar 14% variasi dapat dijelaskan model.
- McKelvey & Zavoina = 0.203 → model cukup baik dalam latent variable explanation.
- Tjur R² = 0.089 → model memisahkan rata-rata probabilitas antar kelas sekitar 9%.
- Secara keseluruhan: model memiliki goodness of fit yang cukup, meskipun belum terlalu kuat.
- McFadden R² = 0.09 → Umumnya > 0.2 baru dianggap “baik”, artinya model Anda masih bisa ditingkatkan.
Bab XII Distribusi Multinomial
Distribusi multinomial merupakan perluasan dari distribusi binomial untuk kasus dengan lebih dari dua kategori. Misalkan \(X_1, X_2, ..., X_k\) adalah jumlah kejadian dalam masing-masing dari \(k\) kategori, maka probabilitas gabungan distribusi multinomial diberikan oleh:
\[ P(X_1 = x_1, ..., X_k = x_k) = \frac{n!}{x_1! x_2! \cdots x_k!} p_1^{x_1} p_2^{x_2} \cdots p_k^{x_k} \]
dengan ketentuan:
- \(\sum_{i=1}^{k} x_i = n\)
- \(\sum_{i=1}^{k} p_i = 1\)
Distribusi ini sangat berguna dalam model klasifikasi dengan lebih dari dua kelas.
12.1 Multinomial Logistic Regression
Model regresi logistik multinomial digunakan untuk menganalisis hubungan antara satu variabel respon kategorik (dengan lebih dari dua kategori) dan satu atau lebih variabel prediktor.
Misalkan variabel respon \(Y\) memiliki \(K\) kategori, dan kita memilih kategori \(K\) sebagai kategori acuan (baseline). Maka, model logit untuk kategori \(j\) (dengan \(j = 1, 2, \dots, K-1\)) adalah:
\[ \log\left(\frac{P(Y = j)}{P(Y = K)}\right) = \beta_{j0} + \beta_{j1}x_1 + \dots + \beta_{jp}x_p \]
12.1.1 Model Logit Kategori Acuan (Baseline-category logit model)
Model logit kategori acuan membandingkan setiap kategori \(j\) dengan kategori acuan \(c\). Bentuk umum model ini:
\[ \log\left(\frac{\pi_j}{\pi_c}\right) = \alpha_j + \beta_j x, \quad j = 1, \dots, c-1 \]
Dengan:
- \(\pi_j\): probabilitas respon berada di kategori \(j\)
- \(\pi_c\): probabilitas respon berada di kategori acuan (baseline)
Maka, jumlah fungsi logit yang terbentuk adalah sebanyak \(c - 1\). Kategori baseline secara default
adalah kategori terakhir, tetapi dapat diubah menggunakan fungsi
relevel() di R.
Contoh Kasus: Respon dengan 3 Kategori
Jika variabel respon \(Y \in \{1, 2, 3\}\) dan kategori ke-3 dipilih sebagai baseline, maka:
\[ \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 \]
Jika ingin menghitung logit antara kategori 1 dan 2:
\[ \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) = (\alpha_1 + \beta_1 x) - (\alpha_2 + \beta_2 x) = (\alpha_1 - \alpha_2) + (\beta_1 - \beta_2)x \]
Ciri-ciri model baseline-category logit:
- Digunakan untuk data kategorik dengan jumlah kategori > 2
- Menghasilkan \((c - 1)\) fungsi logit terhadap satu baseline
- Logit antar kategori lainnya dapat dihitung dari selisih logit terhadap baseline
Implementasi di R dapat menggunakan fungsi multinom()
dari package nnet.
12.1.2 Estimasi Parameter
Parameter dalam model ini diestimasi menggunakan metode maximum likelihood yang biasanya diselesaikan dengan algoritma iteratif seperti Newton-Raphson.
Fungsi Log-Likelihood:
\[ \ell(\beta) = \sum_{i=1}^{n} \sum_{j=1}^{K} y_{ij} \log(\pi_{ij}) \]
Dengan:
- \(\pi_{ij} = P(Y_i = j | x_i)\): probabilitas bahwa observasi ke-\(i\) berada pada kategori \(j\)
- \(y_{ij} = 1\) jika \(Y_i = j\), dan \(0\) jika tidak
12.2 Contoh Kasus
Pemilihan Moda Transportasi Mahasiswa Antar Fakultas
Sebuah universitas di kota besar tengah melakukan studi tentang pemilihan moda transportasi mahasiswa dalam rangka mengembangkan kebijakan transportasi ramah lingkungan dan meningkatkan aksesibilitas ke kampus.
Penelitian ini melibatkan 150 mahasiswa dari tiga fakultas berbeda, yaitu Fakultas Ekonomi, Fakultas Teknik, dan Fakultas Sastra.
Setiap mahasiswa ditanya mengenai moda transportasi utama yang mereka gunakan untuk datang ke kampus: Jalan kaki, Motor, atau Mobil; Jarak tempat tinggal ke kampus (dalam kilometer); Pengalaman kerja paruh waktu selama kuliah (dalam tahun); dan Fakultas tempat mereka belajar.
12.3 Simulasi Data
set.seed(123)
n <- 150
Fakultas <- sample(c("Ekonomi", "Teknik", "Sastra"), n, replace = TRUE)
Jarak <- round(runif(n, 1, 20), 1) # jarak ke kampus dalam km
Experience <- round(pmax(rnorm(n, mean = 7, sd = 3), 0)) # misalnya pengalaman kerja paruh waktu
# Moda transportasi dipilih berdasarkan probabilitas fakultas
Transportasi <- sapply(Fakultas, function(fak) {
if (fak == "Ekonomi") {
sample(c("Motor", "Mobil", "Jalan Kaki"), 1, prob = c(0.6, 0.3, 0.1))
} else if (fak == "Teknik") {
sample(c("Motor", "Mobil", "Jalan Kaki"), 1, prob = c(0.5, 0.4, 0.1))
} else {
sample(c("Motor", "Mobil", "Jalan Kaki"), 1, prob = c(0.3, 0.2, 0.5))
}
})
df <- data.frame(Transportasi = factor(Transportasi), Jarak, Fakultas = factor(Fakultas), Experience)
df$Transportasi <- relevel(df$Transportasi, ref = "Motor")
head(df)12.4 Estimasi Model
library(nnet)
model_mnlogit <- multinom(Transportasi ~ Jarak + Fakultas + Experience, data = df)## # weights: 18 (10 variable)
## initial value 164.791843
## iter 10 value 140.082078
## final value 138.173757
## convergedsummary(model_mnlogit)## Call:
## multinom(formula = Transportasi ~ Jarak + Fakultas + Experience,
## data = df)
##
## Coefficients:
## (Intercept) Jarak FakultasSastra FakultasTeknik Experience
## Jalan Kaki -4.077909 0.04165933 3.5303807 1.6687723 0.03989528
## Mobil -1.180318 0.01156066 0.5148908 0.5233882 0.03979512
##
## Std. Errors:
## (Intercept) Jarak FakultasSastra FakultasTeknik Experience
## Jalan Kaki 1.2898213 0.04407481 1.0674539 1.1305098 0.08246141
## Mobil 0.6624886 0.03390983 0.4903787 0.4477954 0.05879029
##
## Residual Deviance: 276.3475
## AIC: 296.347512.5 Nilai P-Value dan Interpretasi
z <- summary(model_mnlogit)$coefficients / summary(model_mnlogit)$standard.errors
pval <- 2 * (1 - pnorm(abs(z)))
round(pval, 4)## (Intercept) Jarak FakultasSastra FakultasTeknik Experience
## Jalan Kaki 0.0016 0.3446 0.0009 0.1399 0.6285
## Mobil 0.0748 0.7332 0.2937 0.2425 0.4985Interpretasi:
Berdasarkan hasil analisis, terdapat perbedaan pengaruh prediktor pada pilihan moda transportasi:
Untuk pemilihan jalan kaki, variabel Fakultas Sastra berpengaruh signifikan. Mahasiswa dari Fakultas Sastra cenderung lebih memilih jalan kaki dibandingkan mahasiswa dari fakultas referensi.
Untuk pemilihan mobil, tidak ada variabel yang menunjukkan pengaruh signifikan, sehingga pemilihan moda ini kemungkinan dipengaruhi oleh faktor lain yang belum dimasukkan ke dalam model (misalnya, status ekonomi, kepemilikan kendaraan, dsb.).
Jarak tempuh dan Experience tidak menunjukkan pengaruh signifikan dalam pemilihan jalan kaki maupun mobil.
12.6 Prediksi dan Validasi
df$Predicted <- predict(model_mnlogit)
table(Predicted = df$Predicted, Actual = df$Transportasi)## Actual
## Predicted Motor Jalan Kaki Mobil
## Motor 60 12 39
## Jalan Kaki 13 16 9
## Mobil 1 0 012.7 Visualisasi
library(ggplot2)
ggplot(df, aes(x = Jarak, y = Experience, color = Predicted)) +
geom_point(size = 2) +
labs(
title = "Prediksi Moda Transportasi Berdasarkan Jarak dan Pengalaman",
x = "Jarak ke Kampus (km)",
y = "Pengalaman Paruh Waktu (tahun)",
color = "Prediksi Moda"
) +
theme_minimal()
12.8 Kesimpulan
Model regresi logistik multinomial berhasil mengestimasi pengaruh Jarak, Fakultas, dan Experience terhadap pemilihan moda transportasi mahasiswa. Hasil menunjukkan bahwa:
Fakultas mempengaruhi probabilitas memilih transportasi tertentu (misalnya, mahasiswa Sastra cenderung lebih banyak berjalan kaki).
Jarak ke kampus juga berkontribusi: semakin jauh, semakin besar kemungkinan memilih kendaraan bermotor.
Model dapat digunakan untuk prediksi, meskipun evaluasi lanjutan (misalnya cross-validation) diperlukan untuk validasi lebih kuat.
Bab XIII Regresi Logistik Ordinal
Regresi logistik ordinal digunakan ketika variabel respon \(Y\) bersifat ordinal (memiliki urutan), seperti tingkat kepuasan: Rendah, Sedang, Tinggi.
Model ini berbeda dengan:
- Regresi logistik biner: hanya 2 kategori.
- Regresi logistik multinomial: kategori lebih dari 2, tetapi tidak memiliki urutan.
13.1 Model Cumulative Logit
Model yang digunakan adalah Cumulative Logit Model dengan asumsi proportional odds, yaitu:
\[ \log \left( \frac{P(Y \leq j)}{P(Y > j)} \right) = \alpha_j + \beta x \]
- \(\alpha_j\): intercept untuk kategori ke-\(j\).
- \(\beta\): koefisien regresi yang sama untuk semua kategori kumulatif.
Untuk \(c\) kategori, akan ada \(c - 1\) model logit kumulatif.
13.2 Interpretasi Koefisien
Koefisien \(\beta\) menjelaskan pengaruh variabel \(x\) terhadap kemungkinan respon berada pada kategori yang lebih rendah atau sama.
- Jika \(\beta > 0\): semakin besar \(x\), semakin besar peluang berada di kategori rendah.
- Jika \(\beta < 0\): semakin besar \(x\), semakin besar peluang berada di kategori tinggi.
Odds ratio:
\[ \text{OR} = e^{\beta} \]
13.3 Contoh Kasus
Dalam lanjutan studi perangkat yang digunakan oleh mahasiswa dari berbagai jurusan, tim peneliti tertarik melihat apakah usia mahasiswa memengaruhi kepuasan mereka terhadap perangkat (dinyatakan berdasarkan pengalaman kerja sebagai indikator tidak langsung).
# Transformasi Experience menjadi satisfaction ordinal
df$satisfaction <- cut(df$Experience,
breaks = c(-Inf, 4, 8, Inf),
labels = c("Rendah", "Sedang", "Tinggi"),
ordered_result = TRUE)Estimasi Model Ordinal
library(MASS)##
## Attaching package: 'MASS'## The following object is masked from 'package:dplyr':
##
## selectmodel_ord <- polr(satisfaction ~ Jarak, data = df, Hess = TRUE)
summary(model_ord)## Call:
## polr(formula = satisfaction ~ Jarak, data = df, Hess = TRUE)
##
## Coefficients:
## Value Std. Error t value
## Jarak -0.02185 0.02763 -0.7909
##
## Intercepts:
## Value Std. Error t value
## Rendah|Sedang -1.4205 0.3515 -4.0407
## Sedang|Tinggi 0.4388 0.3312 1.3247
##
## Residual Deviance: 320.3063
## AIC: 326.3063Nilai P-Value
(ctable <- coef(summary(model_ord)))## Value Std. Error t value
## Jarak -0.02185288 0.0276313 -0.7908741
## Rendah|Sedang -1.42046555 0.3515365 -4.0407338
## Sedang|Tinggi 0.43876261 0.3312060 1.3247424p <- pnorm(abs(ctable[, "t value"]), lower.tail = FALSE) * 2
(ctable <- cbind(ctable, "p value" = round(p, 4)))## Value Std. Error t value p value
## Jarak -0.02185288 0.0276313 -0.7908741 0.4290
## Rendah|Sedang -1.42046555 0.3515365 -4.0407338 0.0001
## Sedang|Tinggi 0.43876261 0.3312060 1.3247424 0.1853Prediksi Probabilitas
newdata <- data.frame(Jarak = seq(2, 18, by = 4))
predict(model_ord, newdata = newdata, type = "probs")## Rendah Sedang Tinggi
## 1 0.2015299 0.4168007 0.3816694
## 2 0.2159631 0.4227723 0.3612645
## 3 0.2311309 0.4275208 0.3413483
## 4 0.2470282 0.4309952 0.3219766
## 5 0.2636440 0.4331578 0.3031982Interpretasi:
Faktor yang Mempengaruhi Loyalitas: Kepuasan kerja (baik dalam bentuk kategori maupun skor) dan lama bekerja berpengaruh signifikan terhadap loyalitas karyawan. Karyawan dengan tingkat kepuasan kerja yang lebih tinggi dan lama kerja yang lebih panjang cenderung lebih loyal.
Model Terbaik: Model dengan skor kepuasan sebagai prediktor memberikan hasil yang sedikit lebih baik (AIC lebih rendah). Namun, perbedaan performa antar model tidak besar, sehingga model dengan kategori kepuasan tetap relevan secara interpretatif.
Implikasi dan Rekomendasi: Peningkatan kepuasan kerja dapat menjadi strategi utama dalam meningkatkan loyalitas. Perusahaan perlu memperhatikan retensi jangka panjang dan menciptakan lingkungan kerja yang memuaskan untuk mempertahankan karyawan loyal.
13.4 Asumsi Goodness-of-Fit dan Proportional Odds
Model cumulative logit mengasumsikan efek prediktor sama untuk setiap batas kategori. Jika tidak terpenuhi, dapat dipertimbangkan generalized ordinal model.
13.5 Model Alternatif
Jika asumsi proportional odds tidak terpenuhi, alternatif model ordinal lain adalah:
- Adjacent-category logit
- Continuation-ratio (sequential) logit
13.6 Kesimpulan
- Regresi ordinal efektif untuk data respon berurutan.
- Efek prediktor diinterpretasikan sebagai log-odds kumulatif.
- Diimplementasikan di R menggunakan fungsi
polr()dari packageMASS.
Jika dibutuhkan validasi model, gunakan uji devians atau likelihood ratio test.
13.7 Asumsi Paralelisme
Model cumulative logit didasarkan pada asumsi paralelisme (parallel lines assumption), yaitu:
Koefisien regresi \(\beta\) sama untuk semua kategori kumulatif.
Bentuk umum model:
\[ \log \left( \frac{P(Y \leq j)}{P(Y > j)} \right) = \alpha_j + \beta x, \quad j = 1, \ldots, c - 1 \]
- Hanya nilai \(\alpha_j\) yang berbeda.
- \(\beta\) tetap konstan untuk semua kategori.
Visualisasi: semua kurva logit kumulatif memiliki kemiringan yang sama (parallel), hanya posisi (intercept) yang berbeda.
Konsekuensi Jika Asumsi Tidak Terpenuhi
- Efek prediktor berbeda untuk setiap batas kategori.
- Model cumulative logit menjadi tidak valid.
- Gunakan alternatif:
- Generalized Ordinal Logistic Regression
- Partial Proportional Odds Model
Pengujian Asumsi Paralelisme
Dapat diuji dengan:
- Likelihood Ratio Test
- Brant Test (dengan fungsi
brantdi R)
Bab XIV Model Log Linear Dua Arah
Analisis data kategorik memegang peranan penting dalam statistika terapan, karena banyak fenomena di dunia nyata menghasilkan data dalam bentuk kategori. Contohnya meliputi: jenis kelamin, status pekerjaan, tingkat pendidikan, preferensi konsumen, hingga hasil diagnosis medis. Untuk menganalisis data kategori ini, beberapa pendekatan umum yang digunakan antara lain tabel kontingensi, model log-linear, serta model regresi logistik. Masing-masing metode memiliki keunggulan dan keterbatasan, tergantung pada tujuan analisis serta struktur data yang dianalisis.
Tabel Kontingensi
Sebagai langkah awal eksplorasi, tabel kontingensi digunakan untuk memvisualisasikan hubungan antara dua atau lebih variabel kategorik. Misalnya, dalam studi mengenai pengaruh konsumsi obat terhadap kejadian serangan jantung, tabel kontingensi dapat menyajikan jumlah pasien yang mengalami atau tidak mengalami serangan jantung, berdasarkan jenis obat yang dikonsumsi. Tabel ini membantu dalam mengidentifikasi pola awal serta menghitung ukuran asosiasi seperti odds ratio, risk ratio, maupun chi-square statistic untuk menguji independensi antar variabel.
Pengertian Log-Linear Model
Log-linear model adalah suatu model statistik yang digunakan untuk menganalisis hubungan atau asosiasi antara variabel-variabel kategorik dalam tabel kontingensi. Model ini memodelkan logaritma alami dari frekuensi sel dalam tabel sebagai fungsi linier dari efek-efek utama dan interaksi antar variabel.
Secara matematis, model log-linear dapat dituliskan sebagai:
\[ \log(\mu_{ijk\ldots}) = \lambda + \lambda^A_i + \lambda^B_j + \lambda^C_k + \ldots + \lambda^{AB}_{ij} + \lambda^{AC}_{ik} + \lambda^{BC}_{jk} + \ldots \]
di mana:
- \(\mu_{ijk\ldots}\) adalah ekspektasi (mean) frekuensi untuk sel dengan kombinasi kategori \(i\), \(j\), \(k\), dst.
- \(\lambda\) adalah konstanta (intercept).
- \(\lambda^A_i\), \(\lambda^B_j\), dst. adalah efek utama dari masing-masing variabel.
- \(\lambda^{AB}_{ij}\), dst. adalah efek interaksi antar variabel.
Model log-linear termasuk dalam keluarga Generalized Linear Model (GLM) dengan asumsi bahwa frekuensi sel mengikuti distribusi Poisson. Tidak seperti regresi logistik, log-linear model tidak menganggap adanya variabel dependen atau independen, melainkan semua variabel diperlakukan secara simetris. Oleh karena itu, log-linear model sangat cocok digunakan ketika tujuan analisis adalah untuk memahami struktur asosiasi atau independensi antar variabel, bukan untuk melakukan prediksi.
Dalam tabel kontingensi dengan tiga variabel (misalnya: jenis kelamin, status merokok, dan kejadian penyakit paru), log-linear model dapat digunakan untuk mengevaluasi apakah interaksi dua arah antar variabel cukup untuk menjelaskan data, ataukah perlu memasukkan interaksi tiga arah. Pemilihan model dilakukan dengan membandingkan model sederhana dan model kompleks menggunakan Likelihood Ratio Test.
Regresi Logistik
Di sisi lain, regresi logistik adalah metode yang paling umum digunakan ketika terdapat satu variabel kategorik yang secara eksplisit dianggap sebagai variabel dependen (misalnya: kejadian penyakit: ya/tidak). Satu atau lebih variabel kategorik maupun numerik digunakan sebagai prediktor. Model ini memodelkan logit dari probabilitas kejadian (log odds), dan berguna baik dalam studi observasional maupun eksperimental untuk menjelaskan atau memprediksi kemungkinan suatu outcome.
Ekstensi dari regresi logistik meliputi:
- Regresi logistik multinomial → untuk outcome dengan lebih dari dua kategori.
- Regresi logistik ordinal → untuk outcome dengan kategori berurutan.
Ringkasan
Meskipun ketiga pendekatan ini sama-sama digunakan untuk menganalisis data kategorik, fokus analisis dan interpretasinya berbeda:
- Tabel Kontingensi → bersifat deskriptif, menyajikan frekuensi gabungan antar variabel kategorik.
- Model Log-Linear → bersifat eksploratif, memodelkan asosiasi antar variabel kategorik tanpa menetapkan variabel dependen.
- Model Regresi Logistik → bersifat prediktif, memodelkan probabilitas suatu outcome berdasarkan variabel prediktor.
Pemilihan metode yang tepat sangat bergantung pada tujuan analisis:
- Deskripsi awal data → gunakan tabel kontingensi.
- Eksplorasi struktur asosiasi antar variabel → gunakan model log-linear.
- Prediksi outcome kategorik → gunakan model regresi logistik.
Dalam praktik, ketiga pendekatan ini sering digunakan secara komplementer untuk memperoleh pemahaman yang menyeluruh terhadap data kategorik yang dianalisis.
Contoh Kasus:
Dalam studi kesehatan masyarakat, diteliti apakah vaksin baru efektif mencegah infeksi virus. Subjek dibagi menjadi dua kelompok: satu divaksinasi, satu tidak. Hasil dicatat apakah subjek terinfeksi atau tidak.
| Vaksin | Terinfeksi | Tidak Terinfeksi |
|---|---|---|
| Ya | 15 | 85 |
| Tidak | 40 | 60 |
14.1 Tabel Kontingensi, Model Loglinier, Model Regresi Logistik
# Membuat tabel kontingensi
table_data <- matrix(c(15, 85, 40, 60), nrow = 2,
dimnames = list(Vaksin = c("Ya", "Tidak"),
Infeksi = c("Terinfeksi", "Tidak")))
table_data## Infeksi
## Vaksin Terinfeksi Tidak
## Ya 15 40
## Tidak 85 60# Model log-linear
library(MASS)
loglm(~ Vaksin * Infeksi, data = table_data)## Call:
## loglm(formula = ~Vaksin * Infeksi, data = table_data)
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0 0 1
## Pearson 0 0 1# Regresi logistik
data_glm <- data.frame(
Infeksi = c(1, 0, 1, 0),
Vaksin = factor(c("Ya", "Ya", "Tidak", "Tidak")),
Frek = c(15, 85, 40, 60))
model_logit <- glm(Infeksi ~ Vaksin, weights = Frek, family = binomial, data = data_glm)
summary(model_logit)##
## Call:
## glm(formula = Infeksi ~ Vaksin, family = binomial, data = data_glm,
## weights = Frek)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.4055 0.2041 -1.986 0.046993 *
## VaksinYa -1.3291 0.3466 -3.835 0.000125 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 235.27 on 3 degrees of freedom
## Residual deviance: 219.14 on 2 degrees of freedom
## AIC: 223.14
##
## Number of Fisher Scoring iterations: 514.2 Model Saturated
Model Saturated adalah model log-linear yang memasukkan semua efek utama dan semua interaksi yang mungkin di antara variabel-variabel dalam tabel kontingensi. Dengan kata lain, model ini memiliki jumlah parameter yang sama banyaknya dengan jumlah sel yang diamati, sehingga:
- Fitted frequency = Observed frequency untuk setiap sel.
- Tidak ada residual deviance; deviance = 0.
Rumus model saturated:
\[ \log(\mu_{ijk\ldots}) = \lambda + \sum_i \lambda^A_i + \sum_j \lambda^B_j + \sum_k \lambda^C_k + \sum_{i,j} \lambda^{AB}_{ij} + \sum_{i,k} \lambda^{AC}_{ik} + \sum_{j,k} \lambda^{BC}_{jk} + \lambda^{ABC}_{ijk} \]
# Model saturated (dengan interaksi)
model_saturated <- loglm(~ Vaksin * Infeksi, data = table_data)
model_saturated## Call:
## loglm(formula = ~Vaksin * Infeksi, data = table_data)
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0 0 1
## Pearson 0 0 114.3 Model Independent
Model Independent adalah model log-linear yang hanya memasukkan efek utama saja, tanpa memasukkan interaksi antar variabel. Artinya, model ini mengasumsikan bahwa semua variabel saling bebas (independen) satu sama lain.
Rumus model independent:
\[ \log(\mu_{ijk\ldots}) = \lambda + \lambda^A_i + \lambda^B_j + \lambda^C_k + \ldots \]
Jika model independent tidak sesuai dengan data, ini menunjukkan adanya asosiasi atau ketergantungan antar variabel.
# Model independence (tanpa interaksi)
model_independen <- loglm(~ Vaksin + Infeksi, data = table_data)
model_independen## Call:
## loglm(formula = ~Vaksin + Infeksi, data = table_data)
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 16.12336 1 5.934748e-05
## Pearson 15.67398 1 7.525229e-0514.4 Odds Ratio
Odds Ratio (OR) adalah ukuran asosiasi antara dua variabel kategorik. OR menunjukkan berapa kali lebih besar (atau lebih kecil) peluang kejadian di satu grup dibandingkan dengan grup referensi.
Rumus Odds Ratio:
\[ OR = \frac{a / b}{c / d} = \frac{a \times d}{b \times c} \]
Interpretasi:
- OR = 1 → tidak ada asosiasi (odds sama).
- OR > 1 → odds di Group = 1 lebih besar daripada di Group = 0.
- OR < 1 → odds di Group = 1 lebih kecil daripada di Group = 0.
odds_ratio <- (15 / 85) / (40 / 60)
odds_ratio## [1] 0.264705914.5 Estimasi Parameter
Parameter dalam model log-linear diestimasi menggunakan metode Maximum Likelihood Estimation (MLE).
Untuk model berbasis distribusi Poisson, fungsi likelihood adalah:
\[ L(\theta) = \prod_{i,j,k,\ldots} \frac{\mu_{ijk\ldots}^{n_{ijk\ldots}} \exp^{-\mu_{ijk\ldots}}}{n_{ijk\ldots}!} \]
di mana \(n_{ijk\ldots}\) adalah frekuensi observasi pada sel ke-\(ijk\ldots\).
Kemudian log-likelihood akan dimaksimalkan untuk mendapatkan estimasi parameter \(\lambda\). Proses estimasi dilakukan secara iteratif (misalnya, dengan algoritma Iterative Proportional Fitting (IPF) atau Newton-Raphson) hingga diperoleh parameter yang optimal.
logOR <- log((table_data[1,1] * table_data[2,2]) / (table_data[1,2] * table_data[2,1]))
logOR## [1] -1.329136exp(logOR)## [1] 0.264705914.6 Perbandingan Model
# Perbandingan model: deviance
anova(model_independen, model_saturated)## LR tests for hierarchical log-linear models
##
## Model 1:
## ~Vaksin + Infeksi
## Model 2:
## ~Vaksin * Infeksi
##
## Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1 16.12336 1
## Model 2 0.00000 0 16.12336 1 6e-05
## Saturated 0.00000 0 0.00000 0 1e+00Interpretasi:
Berdasarkan uji perbandingan deviance antara model independen dan saturated, diperoleh p-value sebesar 0.00006. Ini menunjukkan bahwa model independen tidak cocok untuk data, dan bahwa terdapat interaksi signifikan antara status vaksinasi dan kejadian infeksi. Artinya, status vaksinasi memengaruhi kemungkinan seseorang terinfeksi virus — vaksinasi memiliki efek yang signifikan secara statistik.
Bab XV Model Log Linear Tiga Arah
Pada pembahasan sebelumnya, kita telah memahami bahwa salah satu tujuan utama dari penggunaan model log-linear adalah untuk mengestimasi parameter-parameter yang menggambarkan hubungan antar variabel kategorik.
Pada bagian ini, kita akan mempelajari model log-linear yang lebih kompleks, yaitu model yang diterapkan pada tabel kontingensi tiga arah. Model ini melibatkan tiga variabel kategorik, sehingga memungkinkan berbagai jenis interaksi, termasuk interaksi tiga arah, yang mencerminkan pengaruh gabungan ketiga variabel secara simultan.
15.1 Model Log-Linear untuk Tabel Tiga Arah
Model log-linear yang mencakup tiga variabel kategorik (misalnya: X, Y, dan Z) dapat dibangun dalam beberapa bentuk, tergantung pada tingkat interaksi yang diinginkan. Berikut adalah berbagai jenis model log-linear tiga arah yang umum digunakan:
Model Saturated
Model saturated memasukkan semua efek utama dan seluruh kemungkinan interaksi (dua arah dan tiga arah):
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} + \lambda^{XYZ}_{ijk} \]
Model ini paling fleksibel, karena memuat interaksi tertinggi (tiga arah), sehingga frekuensi hasil model akan sepenuhnya sesuai dengan frekuensi observasi.
Model Homogeneous Association
Model ini mempertimbangkan semua interaksi dua arah antar variabel, namun tidak memasukkan interaksi tiga arah:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
Model Conditional Association
Model conditional hanya memasukkan sebagian interaksi dua arah, tergantung pada variabel yang dijadikan kondisi.
- Conditional pada X:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \]
- Conditional pada Y:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \]
- Conditional pada Z:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
Model Joint Independence
Model ini digunakan untuk menguji independensi sepasang variabel dari variabel ketiga.
- Independensi antara X & Y:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} \]
- Independensi antara X & Z:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XZ}_{ik} \]
- Independensi antara Y & Z:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{YZ}_{jk} \]
Model Tanpa Interaksi (Mutual Independence)
Model ini hanya memuat efek utama dari ketiga variabel tanpa mempertimbangkan adanya interaksi:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k \]
15.2 Pengujian Interaksi dalam Model Log-Linear Tiga Arah
Dalam analisis model log-linear tiga arah, pengujian interaksi dilakukan untuk menentukan apakah terdapat hubungan antar variabel yang signifikan. Pengujian dilakukan secara bertingkat, mulai dari interaksi paling kompleks hingga yang lebih sederhana.
Tahapan Pengujian:
Pengujian Interaksi Tiga Arah (XYZ):
- Bandingkan model saturated dengan model homogeneous.
- Tujuan: Menilai apakah interaksi tiga arah diperlukan.
Pengujian Interaksi Dua Arah (XY, XZ, YZ):
- Bandingkan model homogeneous dengan model conditional.
- Bandingkan model conditional dengan model joint independence.
- Bandingkan model joint independence dengan model tanpa interaksi.
Setiap tahapan pengujian bertujuan untuk mengevaluasi kecocokan model dan menentukan tingkat interaksi yang paling sesuai dengan struktur data yang diamati.
15.3 Contoh Kasus
Peneliti ingin mengetahui apakah terdapat hubungan antara tingkat kepuasan terhadap layanan, kelompok usia, dan jenis kelamin. Kategori Variabel: - Usia: 1 = Muda, 2 = Dewasa, 3 = Lansia - Jenis Kelamin: 1 = Laki-laki, 2 = Perempuan - Kepuasan: 1 = Puas, 2 = Tidak Puas
15.4 Analisis Log-Linear untuk Tabel Tiga Arah
# Input data sesuai tabel
z.age <- factor(rep(c("1muda", "2dewasa", "3lansia"), each = 4))
x.sex <- factor(rep(c("1M", "2F"), each = 2, times = 3))
y.satis <- factor(rep(c("1puas", "2tdkpuas"), times = 6))
counts <- c(95, 25, 110, 60, 150, 40, 135, 90, 80, 50, 85, 55)
data <- data.frame(
Usia = z.age,
Jenis_Kelamin = x.sex,
Kepuasan = y.satis,
Frekuensi = counts
)
datatable3d <- xtabs(Frekuensi ~ Usia + Jenis_Kelamin + Kepuasan, data = data)
ftable(table3d)## Kepuasan 1puas 2tdkpuas
## Usia Jenis_Kelamin
## 1muda 1M 95 25
## 2F 110 60
## 2dewasa 1M 150 40
## 2F 135 90
## 3lansia 1M 80 50
## 2F 85 55as.array(table3d)## , , Kepuasan = 1puas
##
## Jenis_Kelamin
## Usia 1M 2F
## 1muda 95 110
## 2dewasa 150 135
## 3lansia 80 85
##
## , , Kepuasan = 2tdkpuas
##
## Jenis_Kelamin
## Usia 1M 2F
## 1muda 25 60
## 2dewasa 40 90
## 3lansia 50 55# Model saturated: semua interaksi hingga orde tiga
model_saturated_3d <- loglm(~ Usia * Jenis_Kelamin * Kepuasan, data = table3d)
summary(model_saturated_3d)## Formula:
## ~Usia * Jenis_Kelamin * Kepuasan
## attr(,"variables")
## list(Usia, Jenis_Kelamin, Kepuasan)
## attr(,"factors")
## Usia Jenis_Kelamin Kepuasan Usia:Jenis_Kelamin Usia:Kepuasan
## Usia 1 0 0 1 1
## Jenis_Kelamin 0 1 0 1 0
## Kepuasan 0 0 1 0 1
## Jenis_Kelamin:Kepuasan Usia:Jenis_Kelamin:Kepuasan
## Usia 0 1
## Jenis_Kelamin 1 1
## Kepuasan 1 1
## attr(,"term.labels")
## [1] "Usia" "Jenis_Kelamin"
## [3] "Kepuasan" "Usia:Jenis_Kelamin"
## [5] "Usia:Kepuasan" "Jenis_Kelamin:Kepuasan"
## [7] "Usia:Jenis_Kelamin:Kepuasan"
## attr(,"order")
## [1] 1 1 1 2 2 2 3
## 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 tanpa interaksi tiga arah (interaksi dua arah saja)
model_2way <- loglm(~ Usia * Jenis_Kelamin + Usia * Kepuasan + Jenis_Kelamin * Kepuasan, data = table3d)
summary(model_2way)## Formula:
## ~Usia * Jenis_Kelamin + Usia * Kepuasan + Jenis_Kelamin * Kepuasan
## attr(,"variables")
## list(Usia, Jenis_Kelamin, Kepuasan)
## attr(,"factors")
## Usia Jenis_Kelamin Kepuasan Usia:Jenis_Kelamin Usia:Kepuasan
## Usia 1 0 0 1 1
## Jenis_Kelamin 0 1 0 1 0
## Kepuasan 0 0 1 0 1
## Jenis_Kelamin:Kepuasan
## Usia 0
## Jenis_Kelamin 1
## Kepuasan 1
## attr(,"term.labels")
## [1] "Usia" "Jenis_Kelamin" "Kepuasan"
## [4] "Usia:Jenis_Kelamin" "Usia:Kepuasan" "Jenis_Kelamin:Kepuasan"
## attr(,"order")
## [1] 1 1 1 2 2 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 7.327358 2 0.02563802
## Pearson 7.352629 2 0.02531610# Model independen total (main effects only)
model_indep_3d <- loglm(~ Usia + Jenis_Kelamin + Kepuasan, data = table3d)
summary(model_indep_3d)## Formula:
## ~Usia + Jenis_Kelamin + Kepuasan
## attr(,"variables")
## list(Usia, Jenis_Kelamin, Kepuasan)
## attr(,"factors")
## Usia Jenis_Kelamin Kepuasan
## Usia 1 0 0
## Jenis_Kelamin 0 1 0
## Kepuasan 0 0 1
## attr(,"term.labels")
## [1] "Usia" "Jenis_Kelamin" "Kepuasan"
## attr(,"order")
## [1] 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 34.05363 7 1.683113e-05
## Pearson 32.20068 7 3.727989e-05# Perbandingan model
anova(model_indep_3d, model_2way, model_saturated_3d)## LR tests for hierarchical log-linear models
##
## Model 1:
## ~Usia + Jenis_Kelamin + Kepuasan
## Model 2:
## ~Usia * Jenis_Kelamin + Usia * Kepuasan + Jenis_Kelamin * Kepuasan
## Model 3:
## ~Usia * Jenis_Kelamin * Kepuasan
##
## Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1 34.053627 7
## Model 2 7.327358 2 26.726270 5 0.00006
## Model 3 0.000000 0 7.327358 2 0.02564
## Saturated 0.000000 0 0.000000 0 1.00000Interpretasi:
Analisis model log-linear tiga arah dengan variabel Usia, Jenis Kelamin, dan Kepuasan menunjukkan bahwa penambahan interaksi antar variabel secara bertahap memberikan peningkatan kecocokan model terhadap data.
Pada Model 1, yang hanya memuat efek utama, deviance masih cukup besar (34.05), mengindikasikan bahwa model tersebut kurang sesuai untuk menggambarkan struktur data. Ketika interaksi dua arah ditambahkan (Model 2), terjadi penurunan deviance yang sangat signifikan (∆Deviance = 26.73, p < 0.001), sehingga interaksi dua arah memberikan kontribusi penting dalam menjelaskan hubungan antar variabel.
Penambahan interaksi tiga arah pada Model 3 menghasilkan penurunan deviance lebih lanjut yang juga signifikan (∆Deviance = 7.33, p = 0.026). Ini berarti bahwa interaksi tiga arah antara Usia, Jenis Kelamin, dan Kepuasan memang berperan dalam menjelaskan data, meskipun kontribusi tambahan ini tidak sebesar transisi dari Model 1 ke Model 2.
Dengan demikian, model terbaik yang paling sesuai dengan data adalah Model 3 yang memuat interaksi tiga arah, meskipun jika tujuan analisis lebih menitikberatkan pada model yang lebih sederhana dan mudah ditafsirkan, Model 2 sudah memberikan peningkatan yang sangat baik dan bisa dipertimbangkan. Keputusan akhir mengenai model yang akan dipilih tetap perlu mempertimbangkan konteks penelitian dan kebutuhan interpretasi.
Bab XVI Uji Model Interaksi Tiga Arah (Saturated vs Homogenous)
16.1 Penentuan Kategori Referensi
x.sex <- relevel(x.sex, ref = "2F")
y.satis <- relevel(y.satis, ref = "2tdkpuas")
z.age <- relevel(z.age, ref = "3lansia")16.2 Model Saturated
model_saturated <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age +
x.sex * y.satis * z.age,
family = poisson(link = "log"))
summary(model_saturated)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age + x.sex * y.satis * z.age,
## family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.00733 0.13484 29.719 < 2e-16 ***
## x.sex1M -0.09531 0.19540 -0.488 0.62572
## y.satis1puas 0.43532 0.17305 2.516 0.01188 *
## z.age1muda 0.08701 0.18668 0.466 0.64114
## z.age2dewasa 0.49248 0.17115 2.877 0.00401 **
## x.sex1M:y.satis1puas 0.03469 0.24989 0.139 0.88961
## x.sex1M:z.age1muda -0.78016 0.30797 -2.533 0.01130 *
## x.sex1M:z.age2dewasa -0.71562 0.27257 -2.625 0.00865 **
## y.satis1puas:z.age1muda 0.17082 0.23602 0.724 0.46922
## y.satis1puas:z.age2dewasa -0.02985 0.22015 -0.136 0.89213
## x.sex1M:y.satis1puas:z.age1muda 0.69418 0.37247 1.864 0.06236 .
## x.sex1M:y.satis1puas:z.age2dewasa 0.88161 0.33561 2.627 0.00862 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 1.975e+02 on 11 degrees of freedom
## Residual deviance: -2.931e-14 on 0 degrees of freedom
## AIC: 97.502
##
## Number of Fisher Scoring iterations: 3exp(model_saturated$coefficients)## (Intercept) x.sex1M
## 55.0000000 0.9090909
## y.satis1puas z.age1muda
## 1.5454545 1.0909091
## z.age2dewasa x.sex1M:y.satis1puas
## 1.6363636 1.0352941
## x.sex1M:z.age1muda x.sex1M:z.age2dewasa
## 0.4583333 0.4888889
## y.satis1puas:z.age1muda y.satis1puas:z.age2dewasa
## 1.1862745 0.9705882
## x.sex1M:y.satis1puas:z.age1muda x.sex1M:y.satis1puas:z.age2dewasa
## 2.0020661 2.4147727Interpretasi:
Efek Utama:
Intercept (4.00733): Ini adalah log dari expected count (jumlah kejadian) untuk kategori referensi (misalnya: perempuan, tidak puas, usia referensi — biasanya tua). Jika kita exp(4.00733) ≈ 55, sesuai dengan perhitungan sebelumnya.
x.sex1M (-0.09531): Menjadi laki-laki (dibanding perempuan) mengurangi log-count sebesar 0.095, atau peluang menjadi sekitar 0.91 kali dari baseline. Efek ini lemah (hampir nol, tidak signifikan).
y.satis1puas (0.43532): Bersikap puas meningkatkan log-count sebesar 0.435, atau peluang menjadi 1.55 kali baseline. Efek positif dan cukup nyata.
z.age1muda (0.08701): Usia muda sedikit meningkatkan log-count (+0.087), atau peluang sekitar 1.09 kali baseline. Efek ini lemah.
z.age2dewasa (0.49248): Usia dewasa meningkatkan log-count cukup besar (+0.492), peluang menjadi 1.64 kali baseline. Efek cukup kuat.
Interaksi Dua Arah:
x.sex1M : y.satis1puas (0.03469): Interaksi laki-laki & puas hanya menaikkan log-count sangat kecil (+0.034), peluang sekitar 1.035 kali — hampir tidak ada efek.
x.sex1M : z.age1muda (-0.78016): Interaksi laki-laki & muda mengurangi log-count sebesar -0.780, peluang turun menjadi sekitar 0.46 kali baseline — efek negatif yang kuat.
x.sex1M : z.age2dewasa (-0.71562): Interaksi laki-laki & dewasa juga menurunkan log-count (-0.716), peluang sekitar 0.49 kali baseline — efek negatif.
y.satis1puas : z.age1muda (0.17082): Interaksi puas & muda menaikkan log-count sebesar +0.17, peluang menjadi sekitar 1.18 kali baseline — efek kecil-positif.
y.satis1puas : z.age2dewasa (-0.02985): Interaksi puas & dewasa hampir tidak berdampak (koefisien mendekati nol).
Interaksi Tiga Arah:
x.sex1M : y.satis1puas : z.age1muda (0.69418): Kombinasi laki-laki & puas & muda menaikkan log-count sebesar +0.694, peluang menjadi sekitar 2.00 kali baseline — efek interaksi yang kuat dan positif.
x.sex1M : y.satis1puas : z.age2dewasa (0.88161): Kombinasi laki-laki & puas & dewasa menaikkan log-count sebesar +0.882, peluang sekitar 2.41 kali baseline — efek interaksi sangat kuat.
Goodness-of-Fit
Model saturated yang Anda gunakan sangat cocok dengan data (deviance ≈ 0), karena seluruh interaksi (termasuk 3-way) sudah dimasukkan.
AIC = 97.5 → digunakan untuk membandingkan apakah model ini layak dipilih dibanding model yang lebih sederhana.
16.3 Model Homogenous
Model Homogenous dalam analisis log-linear adalah model yang memuat efek utama (main effects) dan interaksi dua arah (two-way interactions) antar variabel, tanpa memasukkan interaksi tiga arah.
Model ini digunakan untuk menguji apakah hubungan antar pasangan variabel sudah cukup menjelaskan struktur data, tanpa perlu mempertimbangkan efek interaksi yang melibatkan semua variabel secara bersamaan.
Rumus:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
Keterangan:
- \(\mu_{ijk}\) → Ekspektasi (mean) jumlah observasi pada sel ke-(i,j,k).
- \(\lambda\) → Konstanta intercept (efek keseluruhan).
- \(\lambda^X_i\) → Efek utama dari kategori ke-i pada variabel X.
- \(\lambda^Y_j\) → Efek utama dari kategori ke-j pada variabel Y.
- \(\lambda^Z_k\) → Efek utama dari kategori ke-k pada variabel Z.
- \(\lambda^{XY}_{ij}\) → Efek interaksi dua arah antara X dan Y.
- \(\lambda^{XZ}_{ik}\) → Efek interaksi dua arah antara X dan Z.
- \(\lambda^{YZ}_{jk}\) → Efek interaksi dua arah antara Y dan Z.
model_homogenous <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age,
family = poisson(link = "log"))
summary(model_homogenous)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.15446 0.11424 36.365 < 2e-16 ***
## x.sex1M -0.43402 0.15107 -2.873 0.00407 **
## y.satis1puas 0.17998 0.14108 1.276 0.20203
## z.age1muda -0.09254 0.15837 -0.584 0.55897
## z.age2dewasa 0.26675 0.14436 1.848 0.06462 .
## x.sex1M:y.satis1puas 0.58504 0.14122 4.143 3.43e-05 ***
## x.sex1M:z.age1muda -0.33523 0.17275 -1.941 0.05231 .
## x.sex1M:z.age2dewasa -0.14124 0.15860 -0.890 0.37321
## y.satis1puas:z.age1muda 0.47584 0.18165 2.620 0.00880 **
## y.satis1puas:z.age2dewasa 0.35325 0.16536 2.136 0.03266 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 197.4971 on 11 degrees of freedom
## Residual deviance: 7.3274 on 2 degrees of freedom
## AIC: 100.83
##
## Number of Fisher Scoring iterations: 4Interpretasi:
Efek Utama:
Intercept (4.15446): Ini adalah log expected count (log jumlah kejadian) untuk baseline kategori (misal: sex = perempuan, satis = tidak puas, age = kategori referensi, jika ada).
x.sex1M (-0.43402, p=0.00407): Jika jenis kelamin laki-laki dibandingkan perempuan, log count rata-rata menurun sekitar 0.434, dan ini signifikan. Artinya, laki-laki cenderung memiliki jumlah kejadian lebih sedikit dibanding perempuan (dengan asumsi variabel lain tetap).
y.satis1puas (0.17998, p=0.20203): Efek kepuasan puas terhadap jumlah kejadian tidak signifikan secara statistik.
z.age1muda (-0.09254, p=0.559): Usia muda tidak berpengaruh signifikan terhadap jumlah kejadian dibanding baseline.
z.age2dewasa (0.26675, p=0.065): Usia dewasa menunjukkan kecenderungan peningkatan jumlah kejadian (log count naik 0.267), tapi hampir signifikan pada level 0.05 (p ~0.065).
Interaksi Dua Arah:
x.sex1M : y.satis1puas (0.58504, p=3.43e-05): Interaksi yang sangat signifikan antara laki-laki dan puas. Artinya, efek kepuasan “puas” pada jumlah kejadian berbeda tergantung jenis kelamin. Khususnya, laki-laki yang puas memiliki peningkatan log count sekitar 0.585 dibanding baseline.
x.sex1M : z.age1muda (-0.33523, p=0.05231): Ada kecenderungan interaksi negatif antara laki-laki dan usia muda, yang mendekati signifikan.
x.sex1M : z.age2dewasa (-0.33523, p=0.05231): Ada kecenderungan interaksi negatif antara laki-laki dan usia muda, yang mendekati signifikan.
y.satis1puas : z.age1muda (-0.14124, p=0.37321): Interaksi antara puas dan usia dewasa tidak signifikan. Artinya, kombinasi laki-laki dan dewasa tidak berbeda nyata dari baseline.
y.satis1puas : z.age2dewasa (0.35325, p=0.03266): Interaksi antara puas dan usia dewasa juga signifikan positif.
16.4 Langkah-langkah Pengujian
Hipotesis
- \(H_0\): \(\lambda^{XYZ}_{ijk} = 0\) (Tidak ada interaksi tiga arah; model yang terbentuk adalah model homogen)
- \(H_1\): \(\lambda^{XYZ}_{ijk} \neq 0\) (Ada interaksi tiga arah; model yang terbentuk adalah model saturated)
Tingkat Signifikansi
- \(\alpha = 5\%\)
Statistik Uji
# Hitung deviance antara model homogenous dan saturated
Deviance.model <- model_homogenous$deviance - model_saturated$deviance
Deviance.model## [1] 7.327358# Hitung derajat bebas
derajat.bebas <- model_homogenous$df.residual - model_saturated$df.residual
derajat.bebas## [1] 2# Chi kuadrat tabel untuk alpha = 0.05
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel## [1] 5.991465# Keputusan uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan## [1] "Tolak H0"Daerah Penolakan
- Tolak \(H_0\) jika \(\Delta\)Deviance \(> \chi^2_{0.05, df}\)
Keputusan
- Karena \(\Delta\)Deviance \(> \chi^2_{0.05, df}\), maka Tolak \(H_0\)
Interpretasi
- Pada taraf nyata 5%, didapatkan hasil tolak \(H_0\) atau dapat dikatakan bahwa terdapat interaksi tiga arah antara jenis kelamin, kelompok usia, dan tingkat kepuasan.
Bab XVII Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on X)
Model log-linear conditional pada variabel X mempertimbangkan efek utama serta interaksi dua arah antara variabel X dengan Y dan X dengan Z. Namun, model ini tidak mencakup interaksi antara Y dan Z, maupun interaksi tiga arah di antara ketiga variabel tersebut. Pendekatan ini sering digunakan untuk menganalisis hubungan parsial dalam data kategorik, di mana fokus analisis adalah pada keterkaitan dua variabel dengan mengkondisikan pada variabel ketiga.
Secara matematis, model ini dinyatakan sebagai:
\[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \]
Model ini memungkinkan analisis keterkaitan antara X-Y dan X-Z, sambil mengabaikan kemungkinan hubungan langsung antara Y dan Z. Hal ini berguna dalam situasi di mana variabel X dianggap sebagai variabel pengkondisi atau pengontrol.
17.1 Model Homogenous
model_homogenous <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age,
family = poisson(link = "log"))
summary(model_homogenous)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.15446 0.11424 36.365 < 2e-16 ***
## x.sex1M -0.43402 0.15107 -2.873 0.00407 **
## y.satis1puas 0.17998 0.14108 1.276 0.20203
## z.age1muda -0.09254 0.15837 -0.584 0.55897
## z.age2dewasa 0.26675 0.14436 1.848 0.06462 .
## x.sex1M:y.satis1puas 0.58504 0.14122 4.143 3.43e-05 ***
## x.sex1M:z.age1muda -0.33523 0.17275 -1.941 0.05231 .
## x.sex1M:z.age2dewasa -0.14124 0.15860 -0.890 0.37321
## y.satis1puas:z.age1muda 0.47584 0.18165 2.620 0.00880 **
## y.satis1puas:z.age2dewasa 0.35325 0.16536 2.136 0.03266 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 197.4971 on 11 degrees of freedom
## Residual deviance: 7.3274 on 2 degrees of freedom
## AIC: 100.83
##
## Number of Fisher Scoring iterations: 417.2 Model Conditional on X (Jenis Kelamin)
model_conditional_X <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age,
family = poisson(link = "log"))
summary(model_conditional_X)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.98239 0.10076 39.525 < 2e-16 ***
## x.sex1M -0.45669 0.15578 -2.932 0.00337 **
## y.satis1puas 0.47608 0.08893 5.354 8.63e-08 ***
## z.age1muda 0.19416 0.11413 1.701 0.08890 .
## z.age2dewasa 0.47446 0.10764 4.408 1.05e-05 ***
## x.sex1M:y.satis1puas 0.56281 0.14029 4.012 6.03e-05 ***
## x.sex1M:z.age1muda -0.27420 0.17044 -1.609 0.10767
## x.sex1M:z.age2dewasa -0.09497 0.15666 -0.606 0.54438
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 197.497 on 11 degrees of freedom
## Residual deviance: 14.911 on 4 degrees of freedom
## AIC: 104.41
##
## Number of Fisher Scoring iterations: 417.3 Langkah-langkah Pengujian
Hipotesis
- \(H_0\): \(\lambda^{YZ}_{jk} = 0\) (Tidak ada interaksi antara kelompok usia dan tingkat kepuasan terhadap jenis kelamin)
- \(H_1\): \(\lambda^{YZ}_{jk} \neq 0\) (Ada interaksi antara kelompok usia dan tingkat kepuasan terhadap jenis kelamin)
Tingkat Signifikansi
- \(\alpha = 5\%\)
Statistik Uji
# Hitung deviance
Deviance.model <- model_conditional_X$deviance - model_homogenous$deviance
Deviance.model## [1] 7.583399# Hitung derajat bebas
derajat.bebas <- model_conditional_X$df.residual - model_homogenous$df.residual
derajat.bebas## [1] 2# Chi kuadrat tabel untuk alpha = 0.05
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel## [1] 5.991465# Keputusan uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan## [1] "Tolak H0"Daerah Penolakan
- Tolak \(H_0\) jika \(\Delta\)Deviance \(> \chi^2_{0.05, df}\)
Keputusan
- Karena \(\Delta\)Deviance \(> \chi^2_{0.05, df}\), maka Tolak \(H_0\)
Interpretasi
- Pada taraf nyata 5%, didapatkan hasil tolak \(H_0\) atau dapat dikatakan bahwa terdapat interaksi antara kelompok usia dan tingkat kepuasan terhadap jenis kelamin.
Bab XVIII Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on Y)
Model log-linear conditional pada variabel Y menyertakan efek utama dan interaksi dua arah antara Y dengan X, serta Y dengan Z. Dalam model ini, interaksi langsung antara X dan Z tidak dimasukkan, begitu pula interaksi tiga arah. Tujuan dari model ini adalah untuk mengevaluasi hubungan antara X dan Z dengan mengontrol pengaruh dari variabel Y.
Secara umum, bentuk model ini dapat dituliskan sebagai berikut:
\[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \]
Model ini mengasumsikan bahwa hubungan antara X dan Z hanya terjadi secara tidak langsung melalui Y. Ini bermanfaat dalam studi di mana Y merupakan variabel yang diyakini memediasi atau mengkondisikan hubungan antara dua variabel lainnya.
18.1 Model Homogenous
model_homogenous <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age,
family = poisson(link = "log"))
summary(model_homogenous)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.15446 0.11424 36.365 < 2e-16 ***
## x.sex1M -0.43402 0.15107 -2.873 0.00407 **
## y.satis1puas 0.17998 0.14108 1.276 0.20203
## z.age1muda -0.09254 0.15837 -0.584 0.55897
## z.age2dewasa 0.26675 0.14436 1.848 0.06462 .
## x.sex1M:y.satis1puas 0.58504 0.14122 4.143 3.43e-05 ***
## x.sex1M:z.age1muda -0.33523 0.17275 -1.941 0.05231 .
## x.sex1M:z.age2dewasa -0.14124 0.15860 -0.890 0.37321
## y.satis1puas:z.age1muda 0.47584 0.18165 2.620 0.00880 **
## y.satis1puas:z.age2dewasa 0.35325 0.16536 2.136 0.03266 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 197.4971 on 11 degrees of freedom
## Residual deviance: 7.3274 on 2 degrees of freedom
## AIC: 100.83
##
## Number of Fisher Scoring iterations: 418.2 Model Conditional on Y (Kelompok Usia)
model_conditional_Y <- glm(counts ~ x.sex + y.satis + z.age +
y.satis * x.sex + y.satis * z.age,
family = poisson(link = "log"))
summary(model_conditional_Y)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + y.satis * x.sex +
## y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.2086 0.1062 39.632 < 2e-16 ***
## x.sex1M -0.5781 0.1165 -4.962 6.99e-07 ***
## y.satis1puas 0.2118 0.1373 1.543 0.1229
## z.age1muda -0.2113 0.1459 -1.448 0.1475
## z.age2dewasa 0.2136 0.1312 1.628 0.1036
## x.sex1M:y.satis1puas 0.5628 0.1403 4.012 6.03e-05 ***
## y.satis1puas:z.age1muda 0.4284 0.1795 2.386 0.0170 *
## y.satis1puas:z.age2dewasa 0.3330 0.1637 2.034 0.0419 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 197.497 on 11 degrees of freedom
## Residual deviance: 11.168 on 4 degrees of freedom
## AIC: 100.67
##
## Number of Fisher Scoring iterations: 418.3 Langkah-langkah Pengujian
Hipotesis
- \(H_0\): \(\lambda^{XZ}_{ik} = 0\) (Tidak ada interaksi antara jenis kelamin dan tingkat kepuasan terhadap kelompok usia)
- \(H_1\): \(\lambda^{XZ}_{ik} \neq 0\) (Ada interaksi antara jenis kelamin dan tingkat kepuasan terhadap kelompok usia)
Tingkat Signifikansi
- \(\alpha = 5\%\)
Statistik Uji
# Hitung deviance
Deviance.model <- model_conditional_Y$deviance - model_homogenous$deviance
Deviance.model## [1] 3.840192# Hitung derajat bebas
derajat.bebas <- model_conditional_Y$df.residual - model_homogenous$df.residual
derajat.bebas## [1] 2# Chi kuadrat tabel untuk alpha = 0.05
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel## [1] 5.991465# Keputusan uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan## [1] "Terima H0"Daerah Penolakan
- Tolak \(H_0\) jika \(\Delta\)Deviance \(> \chi^2_{0.05, df}\)
Keputusan
- Karena \(\Delta\)Deviance \(< \chi^2_{0.05, df}\), maka Terima \(H_0\)
Interpretasi
- Pada taraf nyata 5%, didapatkan hasil terima \(H_0\) atau dapat dikatakan bahwa tidak ada interaksi antara jenis kelamin dan tingkat kepuasan terhadap kelompok usia.
Bab XIX Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on Z)
Model log-linear conditional pada variabel Z menyertakan efek utama serta interaksi dua arah antara Z dengan X dan Z dengan Y. Dalam model ini, tidak terdapat interaksi langsung antara X dan Y, maupun interaksi tiga arah. Tujuan utamanya adalah mengevaluasi hubungan antara X dan Y ketika dikondisikan terhadap variabel Z.
Secara matematis, model ini dapat dituliskan sebagai:
\[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
Model ini cocok digunakan ketika variabel Z diperlakukan sebagai variabel pengontrol, memungkinkan kita untuk memahami hubungan antara X dan Y setelah mengeliminasi pengaruh dari Z.
19.1 Model Homogenous
model_homogenous <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age,
family = poisson(link = "log"))
summary(model_homogenous)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.15446 0.11424 36.365 < 2e-16 ***
## x.sex1M -0.43402 0.15107 -2.873 0.00407 **
## y.satis1puas 0.17998 0.14108 1.276 0.20203
## z.age1muda -0.09254 0.15837 -0.584 0.55897
## z.age2dewasa 0.26675 0.14436 1.848 0.06462 .
## x.sex1M:y.satis1puas 0.58504 0.14122 4.143 3.43e-05 ***
## x.sex1M:z.age1muda -0.33523 0.17275 -1.941 0.05231 .
## x.sex1M:z.age2dewasa -0.14124 0.15860 -0.890 0.37321
## y.satis1puas:z.age1muda 0.47584 0.18165 2.620 0.00880 **
## y.satis1puas:z.age2dewasa 0.35325 0.16536 2.136 0.03266 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 197.4971 on 11 degrees of freedom
## Residual deviance: 7.3274 on 2 degrees of freedom
## AIC: 100.83
##
## Number of Fisher Scoring iterations: 419.2 Model Conditional on Z (Tingkat Kepuasan)
model_conditional_Z <- glm(counts ~ x.sex + y.satis + z.age +
z.age * x.sex + z.age * y.satis,
family = poisson(link = "log"))
summary(model_conditional_Z)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + z.age * x.sex +
## z.age * y.satis, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.99718 0.11386 35.108 < 2e-16 ***
## x.sex1M -0.07411 0.12180 -0.608 0.542895
## y.satis1puas 0.45199 0.12484 3.621 0.000294 ***
## z.age1muda -0.08861 0.16481 -0.538 0.590807
## z.age2dewasa 0.25818 0.15063 1.714 0.086539 .
## x.sex1M:z.age1muda -0.27420 0.17044 -1.609 0.107673
## x.sex1M:z.age2dewasa -0.09497 0.15666 -0.606 0.544381
## y.satis1puas:z.age1muda 0.42837 0.17952 2.386 0.017022 *
## y.satis1puas:z.age2dewasa 0.33297 0.16366 2.034 0.041903 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 197.497 on 11 degrees of freedom
## Residual deviance: 24.873 on 3 degrees of freedom
## AIC: 116.38
##
## Number of Fisher Scoring iterations: 419.3 Langkah-langkah Pengujian
Hipotesis
- \(H_0\): \(\lambda^{XY}_{ij} = 0\) (Tidak ada interaksi antara jenis kelamin dan kelompok usia terhadap tingkat kepuasan)
- \(H_1\): \(\lambda^{XY}_{ij} \neq 0\) (Ada interaksi antara jenis kelamin dan kelompok usia terhadap tingkat kepuasan)
Tingkat Signifikansi
- \(\alpha = 5\%\)
Statistik Uji
# Hitung deviance
Deviance.model <- model_conditional_Z$deviance - model_homogenous$deviance
Deviance.model## [1] 17.54591# Hitung derajat bebas
derajat.bebas <- model_conditional_Z$df.residual - model_homogenous$df.residual
derajat.bebas## [1] 1# Chi kuadrat tabel untuk alpha = 0.05
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel## [1] 3.841459# Keputusan uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan## [1] "Tolak H0"Daerah Penolakan
- Tolak \(H_0\) jika \(\Delta\)Deviance \(> \chi^2_{0.05, df}\)
Keputusan
- Karena \(\Delta\)Deviance \(> \chi^2_{0.05, df}\), maka Tolak \(H_0\)
Interpretasi
- Pada taraf nyata 5%, didapatkan hasil tolak \(H_0\) atau dapat dikatakan bahwa terdapat interaksi antara jenis kelamin dan kelompok usia terhadap tingkat kepuasan.
Bab XX Pemilihan Model Terbaik
20.1 Ringkasan Model Log Linier
library(knitr)
model_data <- data.frame(
Model = c(
"Model Saturated",
"Model Homogenous",
"Model Conditional on X (Jenis Kelamin)",
"Model Conditional on Y (Kelompok Usia)",
"Model Conditional on Z (Tingkat Kepuasan)"
),
`Null Deviance` = c("197.5", "197.5", "197.5", "197.5", "197.5"),
`Residual Deviance` = c("-2.93e-14", "7.33", "14.91", "11.17", "24.87"),
`df Residual` = c(0, 2, 4, 4, 3),
AIC = c(97.50, 100.83, 104.41, 100.67, 116.38)
)
kable(model_data, caption = "Ringkasan Statistik dari Berbagai Model Log-Linear", align = "lcccc")| Model | Null.Deviance | Residual.Deviance | df.Residual | AIC |
|---|---|---|---|---|
| Model Saturated | 197.5 | -2.93e-14 | 0 | 97.50 |
| Model Homogenous | 197.5 | 7.33 | 2 | 100.83 |
| Model Conditional on X (Jenis Kelamin) | 197.5 | 14.91 | 4 | 104.41 |
| Model Conditional on Y (Kelompok Usia) | 197.5 | 11.17 | 4 | 100.67 |
| Model Conditional on Z (Tingkat Kepuasan) | 197.5 | 24.87 | 3 | 116.38 |
20.2 Ringkasan Pengujian Interaksi
library(knitr)
uji_model <- data.frame(
`Perbandingan Model` = c(
"Saturated vs Homogenous",
"Homogenous vs Conditional on X (Jenis Kelamin)",
"Homogenous vs Conditional on Y (Kelompok Usia)",
"Homogenous vs Conditional on Z (Tingkat Kepuasan)"
),
`Δ Deviance` = c(7.327, 7.583, 3.840, 17.546),
`Δ df` = c(2, 2, 2, 1),
`Chi-kuadrat Tabel (α = 0.05)` = c(5.991, 5.991, 5.991, 3.841),
`Keputusan Uji` = c("Tolak H₀", "Tolak H₀", "Terima H₀", "Tolak H₀")
)
kable(uji_model, caption = "Hasil Uji Chi-Kuadrat Antar Model Log-Linear", align = "lcccc")| Perbandingan.Model | Δ.Deviance | Δ.df | Chi.kuadrat.Tabel..α…0.05. | Keputusan.Uji |
|---|---|---|---|---|
| Saturated vs Homogenous | 7.327 | 2 | 5.991 | Tolak H₀ |
| Homogenous vs Conditional on X (Jenis Kelamin) | 7.583 | 2 | 5.991 | Tolak H₀ |
| Homogenous vs Conditional on Y (Kelompok Usia) | 3.840 | 2 | 5.991 | Terima H₀ |
| Homogenous vs Conditional on Z (Tingkat Kepuasan) | 17.546 | 1 | 3.841 | Tolak H₀ |
20.3 Model Terbaik
Model Conditional on Y (Kelompok Usia) adalah model terbaik karena:
Uji Chi-Kuadrat menunjukkan bahwa tidak ada perbedaan signifikan antara model Homogenous dan Model Conditional on Y (Terima H₀), yang berarti model ini cukup baik menjelaskan data tanpa memerlukan kompleksitas tambahan.
AIC paling rendah di antara model-model yang lebih sederhana (selain model Saturated), yaitu 100.67, sedikit lebih baik dari model Homogenous.
Model Conditional on Y menunjukkan bahwa dengan mengontrol variabel Kelompok Usia, kita bisa cukup baik menjelaskan hubungan antara Jenis Kelamin dan Kepuasan maupun antara Kelompok Usia dan Kepuasan, tanpa perlu menyertakan interaksi langsung antara Jenis Kelamin dan Kepuasan. Ini menjadikannya model yang seimbang antara kesederhanaan dan keakuratan.
Secara umum, bentuk model ini dapat dituliskan sebagai berikut:
\[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \]
20.3.1 Model Conditional on Y (Kelompok Usia)
bestmodel <- glm(counts ~ x.sex + y.satis + z.age +
y.satis * x.sex + y.satis * z.age,
family = poisson(link = "log"))
summary(bestmodel)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + y.satis * x.sex +
## y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.2086 0.1062 39.632 < 2e-16 ***
## x.sex1M -0.5781 0.1165 -4.962 6.99e-07 ***
## y.satis1puas 0.2118 0.1373 1.543 0.1229
## z.age1muda -0.2113 0.1459 -1.448 0.1475
## z.age2dewasa 0.2136 0.1312 1.628 0.1036
## x.sex1M:y.satis1puas 0.5628 0.1403 4.012 6.03e-05 ***
## y.satis1puas:z.age1muda 0.4284 0.1795 2.386 0.0170 *
## y.satis1puas:z.age2dewasa 0.3330 0.1637 2.034 0.0419 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 197.497 on 11 degrees of freedom
## Residual deviance: 11.168 on 4 degrees of freedom
## AIC: 100.67
##
## Number of Fisher Scoring iterations: 4Interpretasi:
- Intercept
Koefisien = 4.2086 → merupakan baseline
log-frekuensi saat: - Jenis kelamin = Perempuan - Kepuasan =
Tidak puas - Usia = Lansia
Dalam kondisi baseline, frekuensi diperkirakan sebesar \(\exp(4.2086) \approx 67.15\).
Nilai ini sangat signifikan (z = 39.632),
artinya intercept berbeda nyata dari nol.
- Jenis Kelamin (Laki-laki) x.sex1M
Koefisien = -0.5781 → responden laki-laki memiliki
log-frekuensi lebih rendah dibanding perempuan.
Rasio frekuensi: \(\exp(-0.5781) \approx
0.561\) → frekuensi laki-laki sekitar 56.1% dari
perempuan.
Signifikan (z = -4.962) → menunjukkan
perbedaan yang bermakna.
- Kepuasan (Puas) y.satis1puas
Koefisien = 0.2118 → responden yang puas memiliki
frekuensi log lebih tinggi dari yang tidak puas.
Rasio: \(\exp(0.2118) \approx 1.236\),
atau peningkatan sekitar 23.6%.
Tidak signifikan (z = 1.543) → belum cukup
bukti bahwa kepuasan berpengaruh langsung.
- Usia Dewasa z.age2dewasa Koefisien =
0.2136 → responden usia dewasa memiliki frekuensi log
lebih tinggi dari lansia.
Rasio: \(\exp(0.2136) \approx 1.238\),
atau kenaikan sekitar 23.8%.
Tidak signifikan (z = 1.628).
- Interaksi: Laki-laki & Puas
x.sex1M:y.satis1puas
Koefisien = 0.5628 → menunjukkan bahwa efek gabungan
laki-laki yang puas menyebabkan peningkatan log-frekuensi.
Rasio: \(\exp(0.5628) \approx 1.756\),
atau naik sekitar 75.6% dibanding model aditif.
Signifikan (z = 4.012) → menunjukkan adanya
interaksi penting antara jenis kelamin dan
kepuasan.
- Interaksi: Puas & Usia Muda
y.satis1puas:z.age1muda Koefisien = 0.4284
→ kombinasi puas dan usia muda meningkatkan log-frekuensi.
Rasio: \(\exp(0.4284) \approx 1.534\),
atau peningkatan sekitar 53.4%.
Signifikan (z = 2.386) → ada hubungan
interaktif signifikan.
- Interaksi: Puas & Usia Dewasa
y.satis1puas:z.age2dewasa Koefisien =
0.3330 → kombinasi puas dan usia dewasa meningkatkan
frekuensi.
Rasio: \(\exp(0.3330) \approx 1.395\),
atau naik sekitar 39.5%.
Signifikan (z = 2.034) → menunjukkan bahwa
efek gabungan signifikan.
Kesimpulan:
- Jenis kelamin (laki-laki) secara signifikan memengaruhi frekuensi.
- Kepuasan dan usia tidak berpengaruh signifikan secara individu.
- Namun, interaksi antar variabel (terutama dengan kepuasan) sangat penting dan signifikan, menunjukkan bahwa hubungan antar variabel lebih kompleks dan tidak bersifat aditif.
20.3.2 Nilai Dugaan Model Terbaik
data.frame(
Usia = z.age,
Jenis_Kelamin = x.sex,
Kepuasan = y.satis,
Frekuensi = counts,
Fitted = bestmodel$fitted.values
)Studi Kasus 1
Pemodelan Tingkat Pengangguran Terbuka di Indonesia menurut Provinsi Tahun 2024 dengan Menggunakan Generalized Linear Model
Latar Belakang
Pengangguran terbuka merupakan salah satu masalah utama yang dihadapi oleh Indonesia sebagai negara berkembang dengan jumlah penduduk yang besar. Pada tahun 2024, Badan Pusat Statistik (BPS) mencatat tingkat pengangguran terbuka (TPT) sebesar 4,91 persen, yang menunjukkan bahwa masih ada sebagian besar angkatan kerja yang belum terserap secara optimal di pasar tenaga kerja[2]. Pengangguran ini berdampak langsung pada kesejahteraan masyarakat dan stabilitas sosial-ekonomi, sehingga menjadi perhatian penting dalam perencanaan pembangunan nasional. Variasi tingkat pengangguran antar provinsi juga menunjukkan adanya ketimpangan yang perlu dianalisis lebih lanjut untuk merumuskan kebijakan yang tepat sasaran.
Faktor-faktor yang mempengaruhi tingkat pengangguran terbuka sangat beragam, mulai dari aspek sosial-kependudukan, pendidikan, hingga kondisi sektor ketenagakerjaan. Data BPS menunjukkan bahwa lulusan sekolah menengah atas dan kejuruan mendominasi kelompok pengangguran, yang disebabkan oleh ketidaksesuaian keterampilan dengan kebutuhan industri dan keterbatasan lapangan kerja[8][4]. Selain itu, persentase penduduk miskin dan pertumbuhan ekonomi yang belum merata juga berkontribusi terhadap tingkat pengangguran di berbagai wilayah[3][6]. Oleh karena itu, pemodelan yang mampu menangkap hubungan kompleks antar variabel tersebut sangat diperlukan untuk memahami faktor-faktor signifikan yang mempengaruhi pengangguran terbuka.
Generalized Linear Model (GLM) menjadi metode yang tepat untuk memodelkan tingkat pengangguran terbuka karena mampu mengakomodasi data dengan distribusi non-normal dan variabel respon kontinu positif. Pendekatan ini memungkinkan analisis yang lebih fleksibel dibanding regresi linier klasik, terutama dalam menangani variabilitas data pengangguran antar provinsi[1]. Dengan menggunakan GLM, penelitian ini bertujuan membangun model yang dapat mengidentifikasi faktor-faktor utama yang mempengaruhi pengangguran terbuka di Indonesia tahun 2024, sekaligus memberikan dasar bagi pengambilan kebijakan yang efektif dalam mengurangi pengangguran secara regional.
Rumusan Masalah
Faktor apa saja yang mempengaruhi tingkat pengangguran terbuka di provinsi-provinsi Indonesia pada tahun 2024?
Bagaimana membangun model prediksi tingkat pengangguran terbuka menggunakan GLM dengan distribusi gamma?
Model sektor mana yang memberikan hasil terbaik berdasarkan kriteria AIC dalam memodelkan pengangguran terbuka?
Tujuan Penelitian
Membentuk model GLM untuk memprediksi tingkat pengangguran terbuka berdasarkan sektor ketenagakerjaan, sosial-kependudukan, dan pendidikan.
Menentukan model terbaik berdasarkan nilai AIC terkecil.
Mengidentifikasi faktor signifikan yang berpengaruh terhadap tingkat pengangguran terbuka di tingkat provinsi.
Data dan Sumber Data
Data yang digunakan oleh penelitian ini bersumber dari Badan Pusat Statistik Tahun 2024. Variabel dalam penelitian ini terdiri atas satu variabel respon dan sembilan variabel prediktor. Pemilihan variabel dilakukan berdasarkan kajian literatur dan ketersediaan data sekunder yang relevan, dengan tujuan untuk menganalisis faktor-faktor yang memengaruhi tingkat pengangguran terbuka.
Variabel respon yang digunakan adalah tingkat pengangguran terbuka, sedangkan variabel prediktor meliputi berbagai aspek ketenagakerjaan, sosial, kependudukan, dan pendidikan. Berikut adalah daftar variabel yang digunakan dalam penelitian ini:
Variabel Respon
Y : Tingkat Pengangguran Terbuka (%)
Variabel Prediktor
Sektor Ketenagakerjaan
X1 : Tingkat Partisipasi Angkatan Kerja (%)
X2 : Jumlah Pencari Kerja Terdaftar (orang)
X3 : Jumlah Lowongan Kerja Terdaftar (lowongan)
Sektor Sosial dan Kependudukan
X4 : Persentase Penduduk Miskin (%)
X5 : Angka Harapan Hidup (tahun)
X6 : Laju Pertumbuhan Penduduk (%)
X7 : Kepadatan Penduduk (jiwa/km²)
Sektor Pendidikan
X8 : Rata-rata Lama Sekolah (tahun)
X9 : Harapan Lama Sekolah (tahun)
Hasil dan Pembahasan
Dengan menggunakan bantuan software Rstudio didapatkan hasil sebagai berikut:
1. Input Data
library(readr)
library(dplyr)
library(broom)
library(MASS)
# Input Data
data <- read.csv("C:/Users/user/Downloads/data adk studi kasus 1.csv", sep=";", dec=",")
head(data)data <- na.omit(data)
# Menyesuaikan nama kolom
colnames(data)[1:11] <- c("Provinsi", "Y", paste0("X", 1:9))
data[ , 1:9] <- lapply(data[ , 1:9], as.numeric)## Warning in lapply(data[, 1:9], as.numeric): NAs introduced by coercion2. Statistika Deskriptif
summary (data)## Provinsi Y X1 X2
## Min. : NA Min. :1.250 Min. :64.63 Min. : 30
## 1st Qu.: NA 1st Qu.:3.206 1st Qu.:67.55 1st Qu.: 2733
## Median : NA Median :4.188 Median :69.85 Median : 8000
## Mean :NaN Mean :4.374 Mean :70.32 Mean : 23939
## 3rd Qu.: NA 3rd Qu.:5.603 3rd Qu.:71.47 3rd Qu.: 16677
## Max. : NA Max. :6.850 Max. :87.88 Max. :296636
## NA's :38
## X3 X4 X5 X6
## Min. : 191 Min. : 3.900 Min. :64.75 Min. :0.310
## 1st Qu.: 1409 1st Qu.: 5.931 1st Qu.:68.87 1st Qu.:1.170
## Median : 5374 Median : 9.850 Median :70.73 Median :1.365
## Mean : 16597 Mean :10.905 Mean :70.58 Mean :1.296
## 3rd Qu.: 10390 3rd Qu.:13.336 3rd Qu.:72.32 3rd Qu.:1.472
## Max. :148663 Max. :31.315 Max. :75.53 Max. :1.930
##
## X7 X8 X9
## Min. : 5.0 Min. : 4.210 Min. : 9.63
## 1st Qu.: 39.5 1st Qu.: 8.300 1st Qu.:12.86
## Median : 100.0 Median : 8.885 Median :13.25
## Mean : 678.2 Mean : 8.841 Mean :13.21
## 3rd Qu.: 251.8 3rd Qu.: 9.515 3rd Qu.:13.72
## Max. :16165.0 Max. :11.490 Max. :15.70
## 3. Uji Signifikansi Koefisien dan Konstanta Berdasarkan Sektor
# Sektor Ketenagakerjaan
model_ket <- glm(Y ~ X1 + X2 + X3, data = data, family = Gamma(link = "log"))
summary(model_ket)##
## Call:
## glm(formula = Y ~ X1 + X2 + X3, family = Gamma(link = "log"),
## data = data)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.095e+00 5.514e-01 11.052 8.47e-13 ***
## X1 -6.677e-02 7.811e-03 -8.548 5.51e-10 ***
## X2 1.642e-07 2.039e-06 0.081 0.936
## X3 2.030e-06 3.398e-06 0.597 0.554
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Gamma family taken to be 0.04335937)
##
## Null deviance: 4.3806 on 37 degrees of freedom
## Residual deviance: 1.4245 on 34 degrees of freedom
## AIC: 101.06
##
## Number of Fisher Scoring iterations: 4AIC(model_ket)## [1] 101.0631# Sektor Sosial dan Kependudukan
model_sos <- glm(Y ~ X4 + X5 + X6 + X7, data = data, family = Gamma(link = "log"))
summary(model_sos)##
## Call:
## glm(formula = Y ~ X4 + X5 + X6 + X7, family = Gamma(link = "log"),
## data = data)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.272e-01 2.464e+00 0.133 0.895
## X4 -1.887e-02 1.165e-02 -1.620 0.115
## X5 1.369e-02 3.181e-02 0.430 0.670
## X6 2.752e-01 2.326e-01 1.183 0.245
## X7 3.099e-05 2.539e-05 1.220 0.231
##
## (Dispersion parameter for Gamma family taken to be 0.1061146)
##
## Null deviance: 4.3806 on 37 degrees of freedom
## Residual deviance: 3.5865 on 33 degrees of freedom
## AIC: 138.51
##
## Number of Fisher Scoring iterations: 6AIC(model_sos)## [1] 138.5117# Sektor Pendidikan
model_pend <- glm(Y ~ X8 + X9, data = data, family = Gamma(link = "log"))
summary(model_pend)##
## Call:
## glm(formula = Y ~ X8 + X9, family = Gamma(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.04211 0.59051 0.071 0.943558
## X8 0.21422 0.05380 3.982 0.000329 ***
## X9 -0.03657 0.06465 -0.566 0.575239
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Gamma family taken to be 0.07079691)
##
## Null deviance: 4.3806 on 37 degrees of freedom
## Residual deviance: 2.6627 on 35 degrees of freedom
## AIC: 123.04
##
## Number of Fisher Scoring iterations: 5AIC(model_pend)## [1] 123.03914. Perbandingan Nilai AIC Ketiga Sektor
aic_values <- data.frame(
Model = c("Ketenagakerjaan", "Sosial-Kependudukan", "Pendidikan"),
AIC = c(AIC(model_ket), AIC(model_sos), AIC(model_pend))
)
aic_values5. Interpretasi Faktor Signifikan
tidy(model_ket) %>% mutate(Sektor = "Ketenagakerjaan") %>%
bind_rows(tidy(model_sos) %>% mutate(Sektor = "Sosial-Kependudukan")) %>%
bind_rows(tidy(model_pend) %>% mutate(Sektor = "Pendidikan")) %>%
filter(p.value < 0.05) %>%
arrange(Sektor)6. Model Final Berdasarkan Variabel Signifikan
# Model ketenagakerjaan
model_ket_sig <- glm(Y ~ X1, data = data, family = Gamma(link = "log"))
summary(model_ket_sig)##
## Call:
## glm(formula = Y ~ X1, family = Gamma(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.197636 0.566828 10.93 5.42e-13 ***
## X1 -0.067661 0.008046 -8.41 5.12e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Gamma family taken to be 0.04610645)
##
## Null deviance: 4.3806 on 37 degrees of freedom
## Residual deviance: 1.6518 on 36 degrees of freedom
## AIC: 102.73
##
## Number of Fisher Scoring iterations: 4# Model sosial dan kependudukan tidak variabel yang signifikan
# Model pendidikan
model_pend_sig <- glm(Y ~ X8, data = data, family = Gamma(link = "log"))
summary(model_pend_sig)##
## Call:
## glm(formula = Y ~ X8, family = Gamma(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.24542 0.31280 -0.785 0.438
## X8 0.19214 0.03505 5.482 3.43e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Gamma family taken to be 0.06852942)
##
## Null deviance: 4.3806 on 37 degrees of freedom
## Residual deviance: 2.6844 on 36 degrees of freedom
## AIC: 121.35
##
## Number of Fisher Scoring iterations: 47. Perbandingan Nilai AIC Model Final
aic_values <- data.frame(
Model = c("Ketenagakerjaan", "Pendidikan"),
AIC = c(AIC(model_ket_sig), AIC(model_pend_sig))
)
aic_values8. Model Terbaik
summary(model_ket_sig)##
## Call:
## glm(formula = Y ~ X1, family = Gamma(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.197636 0.566828 10.93 5.42e-13 ***
## X1 -0.067661 0.008046 -8.41 5.12e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Gamma family taken to be 0.04610645)
##
## Null deviance: 4.3806 on 37 degrees of freedom
## Residual deviance: 1.6518 on 36 degrees of freedom
## AIC: 102.73
##
## Number of Fisher Scoring iterations: 4Sehingga didapatkan model terbaik, yaitu Model Ketenagakerjaan sebagai berikut:
\[ \hat{\mu} = \left( 6.197636 -0.067661 \times X_1 \right)^{-1} \]
9. Cek Asumsi Model Terbaik
model <- model_ket_sig
# 1. Plot Residual Deviance vs Fitted Values
plot(fitted(model), residuals(model, type = "deviance"),
xlab = "Fitted values",
ylab = "Deviance residuals",
main = "Deviance Residuals vs Fitted Values")
abline(h = 0, col = "red", lty = 2)
# 2. Plot Pearson Residuals vs Fitted Values
plot(fitted(model), residuals(model, type = "pearson"),
xlab = "Fitted values",
ylab = "Pearson residuals",
main = "Pearson Residuals vs Fitted Values")
abline(h = 0, col = "blue", lty = 2)
# 3. Cek Overdispersion (tidak umum untuk Gamma tapi boleh dicek)
overdispersion <- sum(residuals(model, type = "pearson")^2) / model$df.residual
print(paste("Overdispersion ratio:", round(overdispersion, 3)))## [1] "Overdispersion ratio: 0.046"# 4. Cek Influential Observations dan Leverage
influence <- influence.measures(model)
summary(influence)## Potentially influential observations of
## glm(formula = Y ~ X1, family = Gamma(link = "log"), data = data) :
##
## dfb.1_ dfb.X1 dffit cov.r cook.d hat
## 38 0.16 -0.16 -0.17 1.95_* 0.01 0.46_*# Visualisasi leverage dan Cook's distance
plot(influence$infmat[, "hat"], ylab = "Leverage (Hat values)", main = "Leverage Plot")
plot(influence$infmat[, "cook.d"], ylab = "Cook's Distance", main = "Cook's Distance Plot")
abline(h = 4 / nrow(model$model), col = "red", lty = 2)
Kesimpulan
Berdasarkan penelitian yang telah dilakukan, didapatkan tiga model yaitu model sektor ketenagakerjaan, model sektor sosial dan kependudukan, serta model sektor pendidikan. Dari masing-masing model tersebut didapatkan model terbaik berdasarkan nilai AIC terkecil yaitu Model Ketenagakerjaan dengan Tingkat Partisipasi Angkatan Kerja sebagai faktor yang memiliki pengaruh secara signifikan terhadap tingkat pengangguran terbuka di Indonesia Tahun 2024.
Selain itu, dapat diketahui bahwa jika variabel lain memiliki nilai konstan, maka tingkat pengangguran terbuka akan berubah dengan sendirinya sebesar nilai konstanta 6,197636. Serta jika variabel lain bernilai konstan, maka Tingkat Partisipasi Angkatan Kerja akan bertambah sebesar -0.067661 untuk setiap satu satuan.
Studi Kasus 2
Pemodelan Tingkat Pengangguran Terbuka di Indonesia menurut Provinsi Tahun 2024 dengan Menggunakan Generalized Linear Model
Latar Belakang
Peran perempuan dalam masyarakat seringkali mengalami konstruksi sosial yang membatasi ruang geraknya, terutama dalam konteks tradisional di mana perempuan dianggap sebagai pengurus rumah tangga dan pendamping suami. Pandangan ini masih melekat kuat dalam berbagai budaya, termasuk di Indonesia, yang berdampak pada ketimpangan peran sosial dan kesempatan perempuan dalam ranah publik. Menurut laporan dari berbagai sumber, termasuk Badan Pusat Statistik dan kajian feminisme, ketidaksetaraan gender tidak hanya disebabkan oleh faktor biologis, tetapi juga oleh struktur sosial dan budaya patriarki yang mengakar kuat dalam masyarakat. Oleh karena itu, memahami sikap masyarakat terhadap peran perempuan sangat penting untuk mendorong perubahan sosial menuju kesetaraan gender.
Survei Social Umum tahun 1975 yang dilakukan oleh National Opinion Research Center (NORC) menyediakan data berharga untuk menganalisis sikap masyarakat terhadap peran wanita, khususnya dalam hal menjaga rumah dan mengikuti suami sepanjang hayat. Data ini mengandung variabel jenis kelamin responden, masa pendidikan, dan jawaban sikap terhadap pernyataan tersebut. Analisis data ini menggunakan tabel kontingensi tiga arah memungkinkan untuk melihat hubungan dan interaksi antar ketiga variabel tersebut secara simultan. Pendekatan ini penting karena sikap sosial tidak hanya dipengaruhi oleh satu faktor saja, melainkan interaksi kompleks antar faktor-faktor sosial yang berbeda.
Penelitian ini akan menggunakan Model Log Linier Tiga Arah untuk menganalisis data tersebut. Model ini memungkinkan pengujian berbagai pola asosiasi antar variabel, seperti model saturated (model lengkap dengan semua interaksi), homogenous association (model dengan interaksi dua arah tanpa interaksi tiga arah), dan conditional association (model dengan interaksi terbatas). Dengan pendekatan ini, penelitian dapat menentukan model yang paling sesuai untuk menggambarkan hubungan antar jenis kelamin, masa pendidikan, dan sikap terhadap peran wanita. Hasilnya diharapkan dapat memberikan pemahaman yang lebih mendalam mengenai faktor-faktor yang mempengaruhi sikap sosial terhadap peran perempuan dalam masyarakat.
Rumusan Masalah
Bagaimana membentuk model log linier tiga arah untuk menganalisis hubungan antar variabel jenis kelamin (A), masa pendidikan (B), dan sikap responden terhadap peran wanita (C)?
Apakah model saturated (model lengkap dengan semua interaksi), homogenous association (model dengan interaksi dua arah tanpa interaksi tiga arah), atau conditional association yang paling cocok untuk data tersebut?
Bagaimana pengujian kecocokan model log linier tiga arah dapat digunakan untuk menentukan model terbaik dalam menjelaskan asosiasi antar variabel?
Tujuan Penelitian
Membentuk dan mengestimasi parameter model log linier tiga arah yang melibatkan variabel jenis kelamin, masa pendidikan, dan sikap responden terhadap peran wanita.
Menguji dan membandingkan kecocokan model saturated, homogenous association, dan conditional association menggunakan statistik uji seperti Chi-square dan nilai AIC untuk memilih model terbaik.
Menjelaskan pola asosiasi antar variabel berdasarkan model log linier yang terpilih untuk memberikan pemahaman yang lebih baik mengenai hubungan antara jenis kelamin, pendidikan, dan sikap sosial.
Data dan Sumber Data
Data yang digunakan dalam penelitian ini berasal dari Pusat Penelitian Pendapat Nasional (National Opinion Research Center - NORC) melalui Survey Social Umum tahun 1975. Survei ini mengumpulkan data mengenai sikap masyarakat terhadap berbagai isu sosial, termasuk sikap terhadap peran wanita dalam rumah tangga. Variabel yang dianalisis meliputi:
A (Jenis Kelamin Responden): pria dan wanita
B (Masa Pendidikan Responden): ≤ 8 tahun, 9-12 tahun, dan ≥ 13 tahun
C (Jawaban Responden): setuju dan tidak setuju terhadap pernyataan “Wanita sebaiknya menjaga atau memelihara rumah mereka dan ikut suami selamanya”
Data disajikan dalam bentuk tabel kontingensi tiga arah yang memungkinkan analisis hubungan antar variabel tersebut.
Hasil dan Pembahasan
Dengan menggunakan bantuan software Rstudio didapatkan hasil sebagai berikut:
1. Input Data
# Input data
counts <- c(
72, 47, # Pria, <=8
110, 196, # Pria, 9-12
44, 179, # Pria, >=13
158, 85, # Wanita, <=8
173, 283, # Wanita, 9-12
28, 187 # Wanita, >=13
)
# Definisi faktor konsisten dengan penamaan sebelumnya
z.age <- factor(rep(c("1<=8", "2:9-12", "3:>=13"), each = 4)) # Pendidikan → z.age
x.sex <- factor(rep(c("1M", "2F"), each = 2, times = 3)) # Jenis Kelamin
y.satis <- factor(rep(c("1puas", "2tdkpuas"), times = 6)) # Kepuasan
# Set kategori referensi
x.sex <- relevel(x.sex, ref = "2F")
y.satis <- relevel(y.satis, ref = "2tdkpuas")
z.age <- relevel(z.age, ref = "3:>=13") # kategori pendidikan tertinggi sebagai referensi
# Gabung jadi data frame
data <- data.frame(
Usia = z.age,
Jenis_Kelamin = x.sex,
Kepuasan = y.satis,
Frekuensi = counts
)
data# Buat tabulasi 3 dimensi
table3d <- xtabs(Frekuensi ~ Usia + Jenis_Kelamin + Kepuasan, data = data)
ftable(table3d)## Kepuasan 2tdkpuas 1puas
## Usia Jenis_Kelamin
## 3:>=13 2F 187 28
## 1M 283 173
## 1<=8 2F 196 110
## 1M 47 72
## 2:9-12 2F 85 158
## 1M 179 442. Model Log-Linear
model_saturated_3d <- loglm(~ Usia * Jenis_Kelamin * Kepuasan, data = table3d)
summary(model_saturated_3d)## Formula:
## ~Usia * Jenis_Kelamin * Kepuasan
## attr(,"variables")
## list(Usia, Jenis_Kelamin, Kepuasan)
## attr(,"factors")
## Usia Jenis_Kelamin Kepuasan Usia:Jenis_Kelamin Usia:Kepuasan
## Usia 1 0 0 1 1
## Jenis_Kelamin 0 1 0 1 0
## Kepuasan 0 0 1 0 1
## Jenis_Kelamin:Kepuasan Usia:Jenis_Kelamin:Kepuasan
## Usia 0 1
## Jenis_Kelamin 1 1
## Kepuasan 1 1
## attr(,"term.labels")
## [1] "Usia" "Jenis_Kelamin"
## [3] "Kepuasan" "Usia:Jenis_Kelamin"
## [5] "Usia:Kepuasan" "Jenis_Kelamin:Kepuasan"
## [7] "Usia:Jenis_Kelamin:Kepuasan"
## attr(,"order")
## [1] 1 1 1 2 2 2 3
## 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 1model_2way <- loglm(~ Usia * Jenis_Kelamin + Usia * Kepuasan + Jenis_Kelamin * Kepuasan, data = table3d)
summary(model_2way)## Formula:
## ~Usia * Jenis_Kelamin + Usia * Kepuasan + Jenis_Kelamin * Kepuasan
## attr(,"variables")
## list(Usia, Jenis_Kelamin, Kepuasan)
## attr(,"factors")
## Usia Jenis_Kelamin Kepuasan Usia:Jenis_Kelamin Usia:Kepuasan
## Usia 1 0 0 1 1
## Jenis_Kelamin 0 1 0 1 0
## Kepuasan 0 0 1 0 1
## Jenis_Kelamin:Kepuasan
## Usia 0
## Jenis_Kelamin 1
## Kepuasan 1
## attr(,"term.labels")
## [1] "Usia" "Jenis_Kelamin" "Kepuasan"
## [4] "Usia:Jenis_Kelamin" "Usia:Kepuasan" "Jenis_Kelamin:Kepuasan"
## attr(,"order")
## [1] 1 1 1 2 2 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 169.9393 2 0
## Pearson 161.4629 2 0model_indep_3d <- loglm(~ Usia + Jenis_Kelamin + Kepuasan, data = table3d)
summary(model_indep_3d)## Formula:
## ~Usia + Jenis_Kelamin + Kepuasan
## attr(,"variables")
## list(Usia, Jenis_Kelamin, Kepuasan)
## attr(,"factors")
## Usia Jenis_Kelamin Kepuasan
## Usia 1 0 0
## Jenis_Kelamin 0 1 0
## Kepuasan 0 0 1
## attr(,"term.labels")
## [1] "Usia" "Jenis_Kelamin" "Kepuasan"
## attr(,"order")
## [1] 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 372.5863 7 0
## Pearson 332.2647 7 0anova(model_indep_3d, model_2way, model_saturated_3d)## LR tests for hierarchical log-linear models
##
## Model 1:
## ~Usia + Jenis_Kelamin + Kepuasan
## Model 2:
## ~Usia * Jenis_Kelamin + Usia * Kepuasan + Jenis_Kelamin * Kepuasan
## Model 3:
## ~Usia * Jenis_Kelamin * Kepuasan
##
## Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1 372.5863 7
## Model 2 169.9393 2 202.6471 5 0
## Model 3 0.0000 0 169.9393 2 0
## Saturated 0.0000 0 0.0000 0 13. Uji Model Interaksi Tiga Arah (Saturated vs Homogenous)
# Model Saturated (semua interaksi 3 arah)
model_saturated <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age +
x.sex * y.satis * z.age,
family = poisson(link = "log"))
summary(model_saturated)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age + x.sex * y.satis * z.age,
## family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.23111 0.07313 71.534 < 2e-16 ***
## x.sex1M 0.41434 0.09424 4.397 1.10e-05 ***
## y.satis1puas -1.89890 0.20264 -9.371 < 2e-16 ***
## z.age1<=8 0.04701 0.10222 0.460 0.6456
## z.age2:9-12 -0.78846 0.13081 -6.027 1.67e-09 ***
## x.sex1M:y.satis1puas 1.40675 0.22445 6.268 3.67e-10 ***
## x.sex1M:z.age1<=8 -1.84231 0.18778 -9.811 < 2e-16 ***
## x.sex1M:z.age2:9-12 0.33040 0.16196 2.040 0.0414 *
## y.satis1puas:z.age1<=8 1.32127 0.23506 5.621 1.90e-08 ***
## y.satis1puas:z.age2:9-12 2.51885 0.24322 10.356 < 2e-16 ***
## x.sex1M:y.satis1puas:z.age1<=8 -0.40260 0.31581 -1.275 0.2024
## x.sex1M:y.satis1puas:z.age2:9-12 -3.42989 0.31110 -11.025 < 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: 5.3722e+02 on 11 degrees of freedom
## Residual deviance: -4.2188e-14 on 0 degrees of freedom
## AIC: 102.02
##
## Number of Fisher Scoring iterations: 3# Model Homogenous (Tanpa interaksi 3 arah)
model_homogeneous <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age,
family = poisson(link = "log"))
summary(model_homogeneous)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.02662 0.07583 66.290 < 2e-16 ***
## x.sex1M 0.73412 0.08907 8.242 < 2e-16 ***
## y.satis1puas -0.89013 0.11387 -7.817 5.41e-15 ***
## z.age1<=8 0.14508 0.10164 1.427 0.153
## z.age2:9-12 -0.08949 0.10468 -0.855 0.393
## x.sex1M:y.satis1puas 0.05966 0.11181 0.534 0.594
## x.sex1M:z.age1<=8 -1.70432 0.13698 -12.442 < 2e-16 ***
## x.sex1M:z.age2:9-12 -0.84589 0.12525 -6.753 1.44e-11 ***
## y.satis1puas:z.age1<=8 0.58432 0.13692 4.268 1.98e-05 ***
## y.satis1puas:z.age2:9-12 0.59384 0.12793 4.642 3.45e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 537.22 on 11 degrees of freedom
## Residual deviance: 169.94 on 2 degrees of freedom
## AIC: 267.96
##
## Number of Fisher Scoring iterations: 5# Hitung deviance dan df
Deviance.model <- model_homogeneous$deviance - model_saturated$deviance
derajat.bebas <- model_homogeneous$df.residual - model_saturated$df.residual
# Chi kuadrat tabel (α = 0.05)
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
# Keputusan uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0 (cukup model 2-way)", "Tolak H0 (butuh interaksi 3-way)")
# Tampilkan hasil
Deviance.model## [1] 169.9393derajat.bebas## [1] 2chi.tabel## [1] 5.991465Keputusan## [1] "Tolak H0 (butuh interaksi 3-way)"4. Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on X)
# Model Homogenous (Tanpa interaksi 3 arah)
model_homogeneous <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age,
family = poisson(link = "log"))
summary(model_homogeneous)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.02662 0.07583 66.290 < 2e-16 ***
## x.sex1M 0.73412 0.08907 8.242 < 2e-16 ***
## y.satis1puas -0.89013 0.11387 -7.817 5.41e-15 ***
## z.age1<=8 0.14508 0.10164 1.427 0.153
## z.age2:9-12 -0.08949 0.10468 -0.855 0.393
## x.sex1M:y.satis1puas 0.05966 0.11181 0.534 0.594
## x.sex1M:z.age1<=8 -1.70432 0.13698 -12.442 < 2e-16 ***
## x.sex1M:z.age2:9-12 -0.84589 0.12525 -6.753 1.44e-11 ***
## y.satis1puas:z.age1<=8 0.58432 0.13692 4.268 1.98e-05 ***
## y.satis1puas:z.age2:9-12 0.59384 0.12793 4.642 3.45e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 537.22 on 11 degrees of freedom
## Residual deviance: 169.94 on 2 degrees of freedom
## AIC: 267.96
##
## Number of Fisher Scoring iterations: 5# Model (kondisional pada x.sex)
model_conditional_X <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age,
family = poisson(link = "log"))
summary(model_conditional_X)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.88054 0.07402 65.935 < 2e-16 ***
## x.sex1M 0.79229 0.09156 8.653 < 2e-16 ***
## y.satis1puas -0.45811 0.07426 -6.169 6.89e-10 ***
## z.age1<=8 0.35295 0.08899 3.966 7.30e-05 ***
## z.age2:9-12 0.12242 0.09363 1.308 0.191
## x.sex1M:y.satis1puas -0.10791 0.10459 -1.032 0.302
## x.sex1M:z.age1<=8 -1.69632 0.13607 -12.466 < 2e-16 ***
## x.sex1M:z.age2:9-12 -0.83774 0.12427 -6.741 1.57e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 537.22 on 11 degrees of freedom
## Residual deviance: 197.72 on 4 degrees of freedom
## AIC: 291.74
##
## Number of Fisher Scoring iterations: 5# Uji Perbandingan Model
# Hitung deviance & derajat bebas
Deviance.model <- model_conditional_X$deviance - model_homogeneous$deviance
derajat.bebas <- model_conditional_X$df.residual - model_homogeneous$df.residual
# Nilai kritis chi-square α = 0.05
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
# Keputusan
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0 (cukup conditional on X)", "Tolak H0 (model conditional on X tidak cukup)")
# Tampilkan hasil
Deviance.model## [1] 27.77594derajat.bebas## [1] 2chi.tabel## [1] 5.991465Keputusan## [1] "Tolak H0 (model conditional on X tidak cukup)"5. Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on Y)
# Model Homogenous (Tanpa interaksi 3 arah)
model_homogeneous <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age,
family = poisson(link = "log"))
summary(model_homogeneous)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.02662 0.07583 66.290 < 2e-16 ***
## x.sex1M 0.73412 0.08907 8.242 < 2e-16 ***
## y.satis1puas -0.89013 0.11387 -7.817 5.41e-15 ***
## z.age1<=8 0.14508 0.10164 1.427 0.153
## z.age2:9-12 -0.08949 0.10468 -0.855 0.393
## x.sex1M:y.satis1puas 0.05966 0.11181 0.534 0.594
## x.sex1M:z.age1<=8 -1.70432 0.13698 -12.442 < 2e-16 ***
## x.sex1M:z.age2:9-12 -0.84589 0.12525 -6.753 1.44e-11 ***
## y.satis1puas:z.age1<=8 0.58432 0.13692 4.268 1.98e-05 ***
## y.satis1puas:z.age2:9-12 0.59384 0.12793 4.642 3.45e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 537.22 on 11 degrees of freedom
## Residual deviance: 169.94 on 2 degrees of freedom
## AIC: 267.96
##
## Number of Fisher Scoring iterations: 5# Model (kondisional pada y.satis)
model_conditional_Y <- glm(counts ~ x.sex + y.satis + z.age +
y.satis * x.sex + y.satis * z.age,
family = poisson(link = "log"))
summary(model_conditional_Y)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + y.satis * x.sex +
## y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.41671 0.05693 95.149 < 2e-16 ***
## x.sex1M 0.08398 0.06404 1.311 0.190
## y.satis1puas -0.79466 0.09942 -7.993 1.32e-15 ***
## z.age1<=8 -0.65967 0.07901 -8.349 < 2e-16 ***
## z.age2:9-12 -0.57678 0.07691 -7.499 6.42e-14 ***
## x.sex1M:y.satis1puas -0.10791 0.10459 -1.032 0.302
## y.satis1puas:z.age1<=8 0.56037 0.12928 4.335 1.46e-05 ***
## y.satis1puas:z.age2:9-12 0.58175 0.12586 4.622 3.80e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 537.22 on 11 degrees of freedom
## Residual deviance: 342.97 on 4 degrees of freedom
## AIC: 436.99
##
## Number of Fisher Scoring iterations: 5# Uji Perbandingan Model
# Hitung deviance & derajat bebas
Deviance.model <- model_conditional_Y$deviance - model_homogeneous$deviance
derajat.bebas <- model_conditional_Y$df.residual - model_homogeneous$df.residual
# Nilai kritis chi-square α = 0.05
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
# Keputusan
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0 (cukup conditional on Y)", "Tolak H0 (model conditional on Y tidak cukup)")
# Tampilkan hasil
Deviance.model## [1] 173.0267derajat.bebas## [1] 2chi.tabel## [1] 5.991465Keputusan## [1] "Tolak H0 (model conditional on Y tidak cukup)"6. Uji Model Interaksi Dua Arah (Homogenous vs Model Conditional on Z)
# Model Homogenous (Tanpa interaksi 3 arah)
model_homogeneous <- glm(counts ~ x.sex + y.satis + z.age +
x.sex * y.satis + x.sex * z.age + y.satis * z.age,
family = poisson(link = "log"))
summary(model_homogeneous)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + x.sex * y.satis +
## x.sex * z.age + y.satis * z.age, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.02662 0.07583 66.290 < 2e-16 ***
## x.sex1M 0.73412 0.08907 8.242 < 2e-16 ***
## y.satis1puas -0.89013 0.11387 -7.817 5.41e-15 ***
## z.age1<=8 0.14508 0.10164 1.427 0.153
## z.age2:9-12 -0.08949 0.10468 -0.855 0.393
## x.sex1M:y.satis1puas 0.05966 0.11181 0.534 0.594
## x.sex1M:z.age1<=8 -1.70432 0.13698 -12.442 < 2e-16 ***
## x.sex1M:z.age2:9-12 -0.84589 0.12525 -6.753 1.44e-11 ***
## y.satis1puas:z.age1<=8 0.58432 0.13692 4.268 1.98e-05 ***
## y.satis1puas:z.age2:9-12 0.59384 0.12793 4.642 3.45e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 537.22 on 11 degrees of freedom
## Residual deviance: 169.94 on 2 degrees of freedom
## AIC: 267.96
##
## Number of Fisher Scoring iterations: 5# Model (kondisional pada y.satis)
model_conditional_Z <- glm(counts ~ x.sex + y.satis + z.age +
z.age * x.sex + z.age * y.satis,
family = poisson(link = "log"))
summary(model_conditional_Z)##
## Call:
## glm(formula = counts ~ x.sex + y.satis + z.age + z.age * x.sex +
## z.age * y.satis, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.01460 0.07272 68.956 < 2e-16 ***
## x.sex1M 0.75185 0.08273 9.088 < 2e-16 ***
## y.satis1puas -0.84943 0.08428 -10.079 < 2e-16 ***
## z.age1<=8 0.14996 0.10158 1.476 0.140
## z.age2:9-12 -0.08978 0.10510 -0.854 0.393
## x.sex1M:z.age1<=8 -1.69632 0.13607 -12.466 < 2e-16 ***
## x.sex1M:z.age2:9-12 -0.83774 0.12427 -6.741 1.57e-11 ***
## y.satis1puas:z.age1<=8 0.56037 0.12928 4.335 1.46e-05 ***
## y.satis1puas:z.age2:9-12 0.58175 0.12586 4.622 3.80e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 537.22 on 11 degrees of freedom
## Residual deviance: 170.22 on 3 degrees of freedom
## AIC: 266.25
##
## Number of Fisher Scoring iterations: 5# Uji Perbandingan Model
# Hitung deviance & derajat bebas
Deviance.model <- model_conditional_Z$deviance - model_homogeneous$deviance
derajat.bebas <- model_conditional_Z$df.residual - model_homogeneous$df.residual
# Nilai kritis chi-square α = 0.05
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
# Keputusan
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0 (cukup conditional on Z)", "Tolak H0 (model conditional on Z tidak cukup)")
# Tampilkan hasil
Deviance.model## [1] 0.2849786derajat.bebas## [1] 1chi.tabel## [1] 3.841459Keputusan## [1] "Terima H0 (cukup conditional on Z)"7. Model Terbaik
aic_saturated <- AIC(model_saturated)
aic_homogeneous <- AIC(model_homogeneous)
aic_conditional_X <- AIC(model_conditional_X)
aic_conditional_Y <- AIC(model_conditional_Y)
aic_conditional_Z <- AIC(model_conditional_Z)
data.frame(
Model = c("Saturated", "Homogeneous", "Conditional on X", "Conditional on Y", "Conditional on Z"),
AIC = c(aic_saturated, aic_homogeneous, aic_conditional_X, aic_conditional_Y, aic_conditional_Z))8. Nilai Dugaan Model Terbaik
data.frame(
Usia = z.age,
Jenis_Kelamin = x.sex,
Kepuasan = y.satis,
Frekuensi = counts,
Fitted = model_conditional_Z$fitted.values
)Kesimpulan
Model Conditional on Z (Jawaban Responden) adalah model terbaik karena:
Uji Chi-Kuadrat menunjukkan bahwa tidak ada perbedaan signifikan antara model Homogenous dan Model Conditional on Z (Terima H₀), yang berarti model ini cukup baik menjelaskan data tanpa memerlukan kompleksitas tambahan.
AIC paling rendah di antara model-model yang lebih sederhana (selain model Saturated), yaitu 266.2474, sedikit lebih baik dari model Homogenous.
Referensi
Greenacre, M. (2007). Correspondence Analysis in Practice. Chapman and Hall/CRC.
McHugh, M. L. (2013). The Chi-square test of independence. Biochemia Medica, 23(2), 143–149.
Agresti, A. (2002). Categorical Data Analysis. Wiley.
Breiman, L. (2001). Random forests. Machine learning, 45(1), 5–32.
Agresti, A. (2018). An Introduction to Categorical Data Analysis (Edisi ke-3). Wiley.
Hosmer, D. W., Lemeshow, S., & Sturdivant, R. X. (2013). Applied Logistic Regression (3rd ed.). Wiley.
Dobson, A. J., & Barnett, A. G. (2008). *An Introduction to Generalized Linear Models*. Chapman & Hall/CRC.
Kleinbaum, D. G., Kupper, L. L., Nizam, A., & Rosenberg, E. S. (2007). Applied Regression Analysis and Other Multivariable Methods. Duxbury Press.
Rothman, K. J., Greenland, S., & Lash, T. L. (2008). Modern Epidemiology. Lippincott Williams & Wilkins.
Hulley, S. B., et al. (2013). Designing Clinical Research (4th ed.). Lippincott Williams & Wilkins.
Newbold, P., Carlson, W., & Thorne, B. (2013). *Statistics for Business and Economics* (8th ed.). Pearson.
McClave, J. T., & Sincich, T. (2014). Statistics (13th ed.). Pearson.
McCullagh, P., & Nelder, J. A. (1989). Generalized Linear Models (2nd ed.). Chapman & Hall.
Breslow, N. E., & Day, N. E. (1980). Statistical Methods in Cancer Research: Volume I. IARC.
Mori, S. A., Martha, S., & Aprizkiyandari, S. (2024). Pemodelan tingkat pengangguran terbuka menggunakan generalized linear model. Buletin Ilmiah Matematika, Statistika dan Terapannya, 13(1), 137–146.
Badan Pusat Statistik Indonesia. (2024). https://www.bps.go.id
Smith, T. W. (1978). General social surveys with comparisons to the omnibus surveys of the Survey Research Center, 1972-1976 (Methodological Report No. 5). National Opinion Research Center