Pendahuluan

Data kategori adalah jenis data yang terdiri dari kelompok atau label tertentu, bukan angka yang bisa diukur secara berkelanjutan. Data ini sering digunakan di berbagai bidang ilmu untuk memahami pola, hubungan, dan tren yang tidak bisa diungkapkan langsung dengan angka. Contoh sederhananya meliputi jenis kelamin (pria/wanita), status pekerjaan (bekerja/tidak bekerja), tingkat pendidikan (SD/SMP/SMA), atau selera pribadi (suka/netral/tidak suka).

Analisis data kategori sangat penting untuk mengeksplorasi dan menafsirkan data yang bersifat nominal (tanpa urutan) atau ordinal (bertingkat). Berbeda dengan data numerik, analisis data kategori membutuhkan metode khusus seperti tabel silang (tabel kontingensi), uji chi-square, regresi logistik, atau model berbasis probabilitas. Seiring perkembangan teknologi dan kecerdasan buatan, teknik analisis data kategori juga semakin banyak dipakai dalam machine learning dan analisis big data.

Tujuan Analisis Data Kategori

Analisis data kategori memiliki beberapa tujuan utama yang mendukung pengambilan keputusan berbasis data. Berikut merupakan tujuan dari analisis data kategori.

Menemukan Pola dan Tren

Analisis data kategori membantu mengungkap pola tersembunyi dalam data. Misalnya, dalam riset pasar, data kategori bisa digunakan untuk melihat preferensi produk berdasarkan usia atau jenis kelamin konsumen.

Menyelidiki Hubungan Antar Variabel

Analisis data kategori membantu mengetahui apakah ada keterkaitan antara dua atau lebih variabel kategori. Contohnya, dalam penelitian kesehatan, analisis ini bisa mengungkap apakah kebiasaan merokok berkaitan dengan risiko penyakit tertentu.

Mendukung Pengambilan Keputusan

Analisis data kategori membantu memahami pola dan hubungan dalam data, sehingga pembuat kebijakan dapat mengambil langkah lebih tepat. Misalnya, pemerintah bisa merancang program bantuan sosial berdasarkan analisis tingkat pendidikan atau pekerjaan masyarakat.

Membangun Model Prediksi

Banyak model statistik menggunakan data kategori sebagai bahan analisis, seperti regresi logistik untuk memprediksi hasil biner. Misalnya, apakah seseorang akan membeli produk atau tidak dalam suatu pusat perbelanjaan.

Definisi dan Ruang Lingkup Analisi Data Kategori

Analisis data kategori adalah metode statistik untuk mengolah variabel yang berbentuk kelompok atau label.

Data Nominal dan Data Ordinal

  • Nominal : Data tanpa urutan (jenis buah, agama).

  • Ordinal : Data dengan urutan (jenjang pendidikan, gaji karyawan).

Data Biner dan Data Multikategori

  • Biner : Memiliki dua kategori (jenis kelamin, status kesehatan).

  • Multikategori : Memiliki lebih dari dua kategori (status sosial, jenjang pendidikan).

Perbedaan dengan Data Kuantitatif

Data kategorik berbeda dengan data kuantitatif karena data kategorik tidak dapat diukur dalam skala numerik kontinu.

Manfaat Analisis Data Kategori dalam Berbagai Bidang

Analisis data kategori memiliki manfaat besar dalam berbagai bidang karena membantu mengungkap pola, tren, dan hubungan antar variabel yang bersifat diskrit atau nominal. Berikut beberapa manfaatnya dalam berbagai bidang :

Bisnis dan Pemasaran

  • Segmentasi Pelanggan : Mengelompokkan pelanggan berdasarkan kategori seperti usia, jenis kelamin, preferensi, atau lokasi.
  • Analisis Preferensi Produk : Melihat produk mana yang paling disukai berdasarkan kategori responden.
  • Pengambilan Keputusan : Menentukan strategi pemasaran berdasarkan perilaku konsumen.

Kesehatan dan Medis

  • Klasifikasi Diagnosis : Mengelompokkan pasien berdasarkan gejala atau hasil tes (positif/negatif).
  • Evaluasi Efektivitas Pengobatan : Menganalisis hubungan antara kategori jenis pengobatan dan hasil pasien.
  • Pemetaan Risiko : Menilai kelompok mana (berdasarkan usia, jenis kelamin, dll.) yang paling berisiko terhadap penyakit tertentu.

Pendidikan

  • Evaluasi Kinerja Mahasiswa : Mengelompokkan berdasarkan nilai mutu (A, B, C, D, E.) atau status kelulusan.

  • Penentuan Jurusan Favorit : Melihat kecenderungan pilihan jurusan berdasarkan kategori wilayah atau gender.

  • Analisis Kepuasan : Mengkategorisasi tanggapan mahasiswa terhadap layanan kampus.

Pemerintah dan Kebijakan Publik

  • Survei Kepuasan Publik : Mengkategorisasi opini masyarakat terhadap layanan publik.
  • Pemilu dan Politik : Menganalisis pemilih berdasarkan kategori wilayah, usia, atau afiliasi partai.
  • Pengambilan Kebijakan : Membuat kebijakan berdasarkan distribusi kategori sosial-ekonomi masyarakat.

Aplikasi dan Media Sosial

  • Analisis Sentimen : Mengkategorisasi komentar atau ulasan (positif, negatif, netral).
  • Targeted Advertising : Mengarahkan iklan ke kategori pengguna tertentu.
  • Pemetaan Perilaku Pengguna : Mengelompokkan pengguna berdasarkan preferensi atau aktivitas.

Teknologi dan Artificial Intelligence

  • Label Encoding : Mengolah data kategori agar bisa digunakan dalam algoritma pembelajaran mesin.

  • Intent Classification : Mengelompokkan pertanyaan pengguna ke dalam kategori permintaan informasi, komplain, permintaan bantuan, dll.

  • Deteksi Intrusi : Mengklasifikasikan aktivitas jaringan sebagai normal atau berbahaya.

Metode dalam Analisis Data Kategori

Terdapat berbagai metode yang dapat digunakan dalam analisis data kategori, tergantung pada tujuan analisis yang ingin dicapai. Berikut merupakan metode-metode yang umum digunakan.

Tabel Kontingensi dan Uji Chi-Square

Metode ini digunakan untuk mengevaluasi apakah terdapat hubungan antara dua variabel kategori. Contohnya adalah menguji keterkaitan antara tingkat pendidikan dengan status pekerjaan.

Regresi Logistik

Metode ini digunakan untuk memprediksi peluang terjadinya suatu peristiwa berdasarkan variabel kategori. Contohnya adalah memprediksi kemungkinan pelanggan melakukan pembelian berdasarkan preferensi kategori produk.

Analisis Correspondence

Metode ini digunakan untuk mengeksplorasi hubungan antar kategori dalam satu dataset. Contohnnya adalah menganalisis preferensi makanan berdasarkan kategori kelompok usia.

Decision Tree dan Random Forest

Metod ini digunakan untuk pembelajaran mesin yang banyak digunakan untuk melakukan klasifikasi berbasis data kategori. Contohnya adalah mengelompokkan pelanggan menurut tingkat loyalitas mereka.

Distribusi Probabilitas dalam Data Kategori

Variabel acak kategori adalah variabel yang hanya dapat memiliki beberapa kategori diskrit sebagai hasilnya. Distribusi probabilitas dari variabel ini menggambarkan kemungkinan terjadinya setiap kategori.

Distribusi Bernoulli

Distribusi Bernoulli digunakan untuk eksperimen biner, yaitu eksperimen yang hanya memiliki dua kemungkinan hasil: sukses (1) dengan peluang sebesar p, dan gagal (0) dengan peluang sebesar 1 - p.

Fungsi probabilitasnya :

\[ P(X = x) = p^x (1 - p)^{1 - x}, \quad x \in \{0, 1\} \]

Penjelasan Notasi :

  • X = Variabel acak biner yang hanya bernilai 0 atau 1.

  • p = Peluang terjadinya sukses (X=1).

Contoh Variabel Acak Bernoulli :

  • Hasil pelemparan koin, misalnya kepala dan ekor.

  • Hasil tendangan pinalti, misalnya gol dan tidak gol.

Perhitungan dengan R :

library(knitr)
## Warning: package 'knitr' was built under R version 4.3.3
library(kableExtra)
## Warning: package 'kableExtra' was built under R version 4.3.3
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.3.3
## 
## Attaching package: 'dplyr'
## The following object is masked from 'package:kableExtra':
## 
##     group_rows
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
set.seed(123)
bernoulli_sample <- rbinom(n = 12, size = 1, prob = 0.5) # 12 percobaan Bernoulli
bernoulli_sample
##  [1] 0 1 0 1 1 0 1 1 1 0 1 0

Distribusi Binomial

Distribusi Binomial merupakan perluasan dari distribusi Bernoulli yang diterapkan pada n kali percobaan yang bersifat independen, di mana setiap percobaan hanya memiliki dua kemungkinan hasil: sukses atau gagal. Bila probabilitas sukses dalam setiap percobaan adalah p.

Fungsi Probabilitasnya :

\[ P(X = k) = \binom{n}{k} p^k (1 - p)^{n - k}, \quad k = 0, 1, 2, ..., n \]

Penjelasan Notasi :

  • X : Jumlah keberhasilan yang diperoleh dari n percobaan.

  • n : Banyaknya percobaan yang dilakukan.

  • k : Banyaknya keberhasilan yang dicatat.

  • p : Peluang sukses dalam satu percobaan.

  • \(\binom{n}{k}\) : Kombinasi “n pilih k”, dihitung dengan rumus \(\dfrac{n!}{k!(n-k)!}\)​

Contoh Variabel Acak Binomial :

  • Jumlah sukses dalam 10 kali pelemparan koin.

  • Jumlah gol dalam 10 tendangan pinalti.

Perhitungan dengan R :

set.seed(123)
binomial_sample <- rbinom(n = 12, size = 6, prob = 0.5) # 12 percobaan Binomial dengan 6 ulangan
binomial_sample
##  [1] 2 4 3 4 5 1 3 5 3 3 5 3

Distribusi Multinomial

Distribusi Multinomial merupakan pengembangan lebih lanjut dari distribusi Binomial, yang digunakan ketika setiap percobaan dapat menghasilkan lebih dari dua hasil yang berbeda. Jika suatu eksperimen dilakukan sebanyak n kali, dan setiap percobaan dapat menghasilkan salah satu dari k kategori dengan peluang masing-masing sebesar p₁, p₂, …, pₖ.

Fungsii Probabilitasnya :

\[ P(X_1 = x_1, ..., X_k = x_k) = \frac{n!}{x_1!x_2!...x_k!} \cdot p_1^{x_1} \cdot p_2^{x_2} \cdots p_k^{x_k} \]

Penjelasan Notasi :

  • Xᵢ : Banyaknya kemunculan kategori ke-i
  • n : Total jumlah percobaan
  • xᵢ : Frekuensi kejadian pada kategori ke-i
  • pᵢ : Probabilitas terjadinya kategori ke-i

Contoh Variabel Acak Multinomial :

  • Distribusi agama siswa dalam satu kelas acak.

  • Distribusi warna permen dalam satu bungkus acak.

Perhitungan dengan R :

set.seed(123)
multinomial_sample <- rmultinom(n = 1, size = 12, prob = c(0.4, 0.5, 0.2))
multinomial_sample
##      [,1]
## [1,]    3
## [2,]    5
## [3,]    4

Distribusi Poisson

Distribusi Poisson digunakan untuk memperkirakan banyaknya kejadian yang terjadi dalam rentang waktu atau ruang tertentu, dengan asumsi rata-rata kejadian sebesar λ per satuan waktu atau ruang.

Fungsi Probabilitasnya :

\[ P(X = k) = \frac{e^{-\lambda} \lambda^k}{k!}, \quad k = 0, 1, 2, ... \]

Penjelasan Notasi :

  • X : Banyaknya kejadian dalam suatu interval tertentu.
  • λ : Rata-rata jumlah kejadian dalam interval tersebut.
  • k : Jumlah kejadian aktual yang diamati.

Contoh Variabel Acak Poisson :

  • Banyaknya panggilan masuk ke pusat layanan pelanggan dalam satu jam.

  • Banyaknya kecelakaan lalu lintas yang terjadi di satu ruas jalan dalam sehari.

Perhitungan dengan R :

set.seed(123)
poisson_sample <- rpois(12, lambda = 4) # 12 sampel dengan rata-rata kejadian 4
poisson_sample
##  [1] 3 6 3 6 7 1 4 7 4 4 8 4

Desain Sampling dalam Analisis Data Kategori

Desain sampling adalah bagian dari tahap awal dalam proses penelitian, khususnya dalam pengumpulan data. Tujuannya adalah agar data yang dikumpulkan bisa digunakan untuk menarik kesimpulan yang akurat tentang populasi.

Dalam analisis data kategori, desain sampling merujuk pada cara atau metode yang digunakan untuk memilih sampel dari populasi agar data yang diperoleh dapat mewakili keseluruhan populasi dengan baik. Ini sangat penting karena pemilihan sampel akan memengaruhi validitas dan reliabilitas hasil analisis.

Secara umum, desain sampling dalam analisis data kategori dapat diklasifikasikan ke dalam dua pendekatan utama, yaitu prospective sampling dan retrosprective sampling. Masing - masing pendekatan ini memiliki karakteristik dan metode sampling yang berbeda.

Prospective Sampling

Prospective sampling adalah metode pengambilan sampel di mana peneliti mengikuti partisipan dari waktu sekarang ke masa depan untuk mengamati kejadian atau hasil tertentu. Ini biasa digunakan dalam penelitian longitudinal atau studi kohort prospektif. Beberapa jenis desain sampling dalam metode ini meliputi :

Eksperimen

Dalam studi eksperimental, subjek secara acak dialokasikan ke dalam kelompok perlakuan dan kontrol. Teknik sampling yang umum digunakan meliputi:

  • Simple Random Sampling : Setiap individu dalam populasi memiliki probabilitas yang sama untuk dipilih.

  • Stratified Random Sampling : Populasi dibagi menjadi strata berdasarkan karakteristik tertentu, lalu sampel diambil secara acak dari setiap strata.

  • Cluster Sampling : Populasi dibagi menjadi kelompok - kelompok (cluster), kemudian beberapa cluster dipilih secara acak untuk dianalisis.

Studi Kohort

Studi kohort adalah jenis penelitian observasional dimana sekelompok orang (kohort) yang memiliki karakteristik tertentu diikuti dalam jangka waktu tertentu untuk melihat bagaimana paparan tertentu memengaruhi kejadian suatu outcome (misalnya penyakit, perilaku, atau kondisi tertentu). Jenis sampling yang umum dalam studi kohort meliputi :

  • Census Sampling : Seluruh anggota dalam populasi tertentu diikutsertakan dalam penelitian.

  • Systematic Sampling : Subjek dipilih berdasarkan interval tertentu dari daftar populasi..

  • Matched Sampling : Setiap individu dalam kelompok kohort dipasangkan dengan individu serupa dalam kelompok lain berdasarkan variabel tertentu.

Retrospective Sampling

Retrospective sampling adalah metode pengambilan sampel dalam penelitian di mana peneliti melihat ke belakang untuk mengidentifikasi subjek berdasarkan outcome yang sudah terjadi, lalu menelusuri paparan atau faktor penyebabnya sebelumnya.

Studi Kasus - Kontrol

Dalam studi kasus-kontrol, sekelompok individu dengan kondisi tertentu dibandingkan dengan kelompok tanpa kondisi tersebut. Teknik sampling yang sering digunakan meliputi :

  • Purposive Sampling : Pemilihan sampel berdasarkan karakteristik yang relevan dengan tujuan penelitian.

  • Snowball Sampling : Responden awal membantu merekrut subjek lain yang memiliki karakteristik serupa.

  • Incidence Density Sampling : Kasus dan kontrol dipilih dari populasi yang sama dengan memperhitungkan periode waktu kemunculan kasus.

Studi Kohort Retrospectif

Dalam studi kohort retrospektif, data historis digunakan untuk mengelompokkan individu berdasarkan paparan dan kemudian menganalisis hasil yang terjadi. Teknik sampling yang sering digunakan meliputi :

  • Convience Sampling : Subek dipilih berdasarkan ketersediaan data yang sudah ada.

  • Quota Sampling : Sampel dipilih untuk mencerminkan proporsi tertentu dalam populasi.

  • Case-Based Sampling : Sampel dipilih berdasarkan karakteristik kasus yang telah terjadi.

Tabel Perbandingan Desain Sampling

Jenis Studi Pendekatan Metode Samping Keuntungan Kelemahan
Eksperimen Prospective SRS, Stratified, Cluster Kontrol tinggi terhadap variabel, hubungan sebab akibat dapat dianalisis Biaya tinggi, etika, dan validitas perlu diperhatikan
Studi Kohort Prospective Census, Systematic, Matched Dapat mengamati perkembangan kejadian dalam jangka panjang Membutuhkan waktu lama, risiko kehilangan partisipan
Studi Kasus - Kontrol Retrospective Purposive, Snowball, Incidnece Density Mudah dan cepat dilakukan, efisien untuk penyakit langka Sulit mengontrol variabel penganggu, rentan bias recall
Studi Kohort Retrospektif Retrospective Convenience, Quota, Case - Based Memanfaatkan data historis, lebih murah daripada studi kohort Ketergantungan pada kualitas data historis, bisa terjadi missing data

Desain pengambilan sampel dalam analisis data kategori sangat dipengaruhi oleh pendekatan yang diterapkan, baik secara prospektif maupun retrospektif. Pemilihan teknik sampling yang tepat, baik dalam eksperimen, studi kohort, maupun studi kasus-kontrol, berperan penting dalam meningkatkan validitas hasil penelitian dan menjamin bahwa temuan dapat diterapkan secara lebih luas. Oleh karena itu, pemahaman yang mendalam mengenai karakteristik tiap metode sampling sangatlah krusial untuk menyusun rancangan penelitian yang kuat dan bermutu tinggi.

Tabel Kontingensi 2 x 2

Tabel kontingensi 2×2 merupakan bentuk paling sederhana dari tabel kontingensi yang digunakan untuk menganalisis hubungan antara dua variabel kategorik. Tabel ini sering digunakan dalam analisis statistik untuk mengevaluasi apakah ada keterkaitan antara dua variabel, seperti jenis pengobatan dengan tingkat kesembuhan atau metode belajar dan hasil kelulusan. Berikut merupakan struktur tabel kontingensi 2 x 2.

Kategori 1 (+) Kategori 2 (-) Total
Grup 1 n11 n12 n1.
Grup 2 n21 n22 n2.
Total n.1 n.2 n

Penjelasan Notasi :

Distribusi Peluang dalam Tabel Kontingensi 2 x 2

Peluang Bersama

Peluang bersama adalah probabilitas bahwa kedua variabel terjadi secara bersamaan dalam suatu sel tabel kontingensi:

\[ P(A_i, B_j) = \frac{n_{ij}}{n} \]

Peluang Marginal

Peluang marginal adalah probabilitas kejadian suatu variabel tanpa mempertimbangkan variabel lainnya.

  • Peluang marginal baris:

\[ P(A_i) = \frac{n_{i.}}{n} \]

  • Peluang marginal kolom:

\[ P(B_j) = \frac{n_{.j}}{n} \]

Peluang Bersyarat

Peluang bersyarat adalah probabilitas suatu kejadian terjadi dengan syarat kejadian lain telah terjadi:

\[ P(B_j | A_i) = \frac{P(A_i, B_j)}{P(A_i)} = \frac{n_{ij}}{n_{i.}} \]

Penerapan dalam Studi Kasus

Tabel berikut merupakan data studi kasus tentang hubungan metode belajar dengan hasil tes siswa di suatu sekolah.

Lulus (+) Tidak Lulus (-) Total
Bimbingan Belajar 80 20 100
Belajar Mandiri 50 50 100
Total 130 70 200

Perhitungan Manual

Langkah 1 : Hitung Peluang Bersama

  • \[P(\text{Bimbingan Belajar, Lulus}) = \frac{80}{200} = 0.4\]
  • \[P(\text{Bimbingan Belajar, Tidak Lulus}) = \frac{20}{200} = 0.1\]
  • \[P(\text{Belajar Mandiri, Lulus}) = \frac{50}{200} = 0.25\]
  • \[P(\text{Belajar Mandiri, Tidak Lulus}) = \frac{50}{200} = 0.25\]

Langkah 2 : Hitung Peluang Marginal

  • \[P(\text{Bimbingan Belajar}) = \frac{100}{200} = 0.5\]
  • \[P(\text{Belajar Mandiri}) = \frac{100}{200} = 0.5\]
  • \[P(\text{Lulus}) = \frac{130}{200} = 0.65\]
  • \[P(\text{Tidak Lulus}) = \frac{70}{200} = 0.35\]

Langkah 3 : Hitung Peluang Bersyarat

  • \[ P(\text{Lulus}|\text{Bimbingan Belajar}) = \frac{80}{100} = 0.8 \]

  • \[ P(\text{Lulus}|\text{Belajar Mandiri}) = \frac{50}{100} = 0.5 \]

  • \[ P(\text{Tidak Lulus}|\text{Bimbingan Belajar}) = \frac{20}{100} = 0.2 \]

  • \[ P(\text{Tidak Lulus}|\text{Belajar Mandiri}) = \frac{50}{100} = 0.5 \]

Perhitungan dengan R

#Data
data <- matrix(c(80,20,50,50), nrow = 2, byrow = TRUE)
colnames(data) <- c("Lulus", "Tidak Lulus")
rownames(data) <- c("Bimbingan Belajar","Belajar Mandiri")
n <- sum(data)

#Peluang Bersama
p_bersama <- data/n

#Peluang Marginal
p_marginal_baris <- rowSums(data)/n
p_marginal_kolom <- colSums(data)/n

#Peluang Bersyarat
p_bersyarat <- data / rowSums(data)

#Hasil
list(Peluang_Bersama = p_bersama, Peluang_Marginal_Baris = p_marginal_baris, Peluang_Marginal_Kolom = p_marginal_kolom, Peluang_Bersyarat = p_bersyarat)
## $Peluang_Bersama
##                   Lulus Tidak Lulus
## Bimbingan Belajar  0.40        0.10
## Belajar Mandiri    0.25        0.25
## 
## $Peluang_Marginal_Baris
## Bimbingan Belajar   Belajar Mandiri 
##               0.5               0.5 
## 
## $Peluang_Marginal_Kolom
##       Lulus Tidak Lulus 
##        0.65        0.35 
## 
## $Peluang_Bersyarat
##                   Lulus Tidak Lulus
## Bimbingan Belajar   0.8         0.2
## Belajar Mandiri     0.5         0.5

Interpretasi :

Peluang Bersama :

  • Peluang siswa mengikuti bimbingan belajar dan lulus adalah 0.40 atau 40%.

  • Peluang siswa mengikuti bimbingan belajar dan tidak lulus adalah 0.10 atau 10%.

  • Peluang siswa belajar mandiri dan lulus adalah 0.25 atau 25%.

  • Peluang siswa belajar mandiri dan tidak lulus adalah 0.25 atau 25%

Peluang Marginal :

  • 50% siswa mengikuti bimbingan belajar, 50% belajar mandiri.
  • 65% siswa lulus, dan 35% tidak lulus.

Peluang Bersyarat :

  • Dari siswa yang mengikuti bimbingan belajar, 80% lulus.
  • Dari siswa yang belajar mandiri, hanya 50% yang lulus.
  • Dapat disimpulkan bahwa bimbingan belajar memiliki tingkat kelulusan yang lebih tinggi dibanding belajar mandiri.

Ukuran Asosiasi dalam Data Kategori 2 x 2

Ukuran asosiasi dalam analisis data kategorik digunakan untuk mengukur kekuatan hubungan antara dua variabel kategorik. Ini penting untuk mengetahui apakah hubungan yang ada signifikan dan seberapa kuat hubungan tersebut. Berikut adalah beberapa ukuran asosiasi yang umum digunakan :

Risk Difference (RD)

Risk Difference (juga dikenal sebagai Absolute Risk Reduction) mengukur selisih risiko antara dua kelompok. Nilai RD menunjukkan perbedaan proporsi kejadian antara kelompok yang terpapar dan tidak terpapar.

Rumus :

\[ RD = \frac{a}{a + b} - \frac{c}{c + d} \]

Dengan :

  • 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 dibandingkan Grup 2.

  • Jika RD=0, maka tidak ada perbedaan risiko antara dua kelompok.

Relative Risk (RR)

Relative Risk membandingkan kemungkinan terjadinya suatu kejadian (misalnya penyakit) pada kelompok yang terpapar terhadap kelompok yang tidak terpapar.

Rumus :

\[ RR = \frac{a / (a + b)}{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 dibandingkan Grup 2.

  • Jika RR=1, maka tidak ada perbedaan risiko antara dua kelompok.

Odds Ratio (OR)

Odds Ratio digunakan untuk membandingkan peluang (odds) terjadinya suatu kejadian antara dua kelompok. Sangat umum dalam studi kasus-kontrol.

Rumus:

\[ 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.

Penerapan dalam Studi Kasus

Tabel berikut merupakan data studi kasus tentang hubungan metode belajar dengan hasil tes siswa di suatu sekolah.

Lulus (+) Tidak Lulus (-) Total
Bimbingan Belajar 80 20 100
Belajar Mandiri 50 50 100
Total 130 70 200

Perhitungan Manual

Langkah 1 : Hitung Risk Difference

  • Risk Difference \[ RD = \frac{a}{a + b} - \frac{c}{c + d} = \frac{80}{100} - \frac{50}{100} = 0.8 - 0.5 = 0.3 \]

Langkah 2 : Hitung Relative Risk

  • Relative Risk \[ RR = \frac{a / (a + b)}{c / (c + d)} = \frac{80/100}{50/100} = 0.8 / 0.5 = 1.6 \]

Langkah 3 : Hitung Odds Ratio

  • Odds Ratio \[ OR = \frac{a \cdot d}{b \cdot c} = \frac{80 \cdot 50}{20 \cdot 50} = \frac{4000}{1000} = 4 \]

Perhitungan dengan R

# Data
a <- 80
b <- 20
c <- 50
d <- 50

# Risk Difference
RD <- (a / (a + b)) - (c / (c + d))

# Relative Risk
RR <- (a / (a + b)) / (c / (c + d))

# Odds Ratio
OR <- (a * d) / (b * c)

#Hasil
list(Risk_Difference = RD, Relative_Risk = RR, Odds_Ratio = OR)
## $Risk_Difference
## [1] 0.3
## 
## $Relative_Risk
## [1] 1.6
## 
## $Odds_Ratio
## [1] 4

Interpretasi :

Risk Difference :

  • Perbedaan risiko (probabilitas) kelulusan antara siswa yang ikut bimbingan belajar dan yang belajar mandiri adalah 0.3 atau 30%. Artinya, siswa bimbingan belajar memiliki tingkat kelulusan 30% lebih tinggi dibanding siswa belajar mandiri.

Relative Risk :

  • Resiko relatif kelulusan antara siswa yang ikut bimbingan belajar dan yang belajar mandiri adalah 1,6 atau 160%. Artinya, siswa bimbingan belajar 1.6 kali lebih mungkin lulus dibandingkan dengan yang belajar mandiri.

Odds Ratio :

  • Odds Ratio kelulusan antara siswa yang ikut bimbingan belajar dan yang belajar mandiri adalah 4 atau 400%. Artinya, peluang relatif untuk lulus 4 kali lebih tinggi pada siswa bimbingan belajar dibandingkan dengan siswa belajar mandiri.

Inferensi Tabel Kontingensi Dua Arah

Inferensi dalam statistik merujuk pada proses menarik kesimpulan tentang suatu populasi berdasarkan data dari sampel. Dalam konteks tabel kontingensi dua arah, inferensi digunakan untuk mengevaluasi hubungan antara dua variabel kategorik yang ditampilkan dalam bentuk tabel. Proses inferensi ini umumnya terbagi menjadi dua jenis utama, yaitu estimasi dan uji hipotesis.

Estimasi

Estimasi bertujuan untuk memperkirakan parameter populasi berdasarkan data sampel. Estimasi dibagi menjadi estimasi titik dan estimasi interval.

Estimasi Titik

Estimasi titik digunakan untuk menentukan satu nilai spesifik sebagai perkiraan terbaik dari parameter populasi.

\[ \hat{p} = \frac{x}{n} \]

  • \(\hat{p}\) adalah estimasi titik proporsi,
  • \(x\) adalah jumlah individu dalam kategori tertentu,
  • \(n\) adalah total jumlah individu dalam sampel.

Estimasi Interval

Estimasi interval bertujuan untuk memberikan rentang nilai yang diyakini mengandung parameter populasi dengan tingkat kepercayaan tertentu.

\[ \hat{p} \pm Z_{\alpha/2} \sqrt{ \frac{\hat{p}(1 - \hat{p})}{n} } \]

  • \(Z_{\alpha/2}\) adalah nilai dari distribusi normal standar untuk tingkat kepercayaan tertentu,
  • \(\hat{p}\) adalah estimasi titik proporsi,
  • \(n\) adalah ukuran sampel.

Uji Hipotesis

Uji Proporsi

Uji proporsi digunakan untuk membandingkan proporsi kejadian antara dua kelompok dalam tabel kontingensi, terutama untuk menentukan apkah terdapat perbedaan yang signifikan dalam proporsi kejadian antara dua kelompok yang berbeda.

Untuk menguji hipotesis bahwa tidak ada perbedaan proporsi antara dua kelompok, kita menggunakan uji z dua proporsi, dengan hipotesis:

- Hipotesis Nol (H0) : tidak ada perbedaan proporsi antara dua kelompok

- Hipotesis Alternatif (H1) : terdapat perbedaan proporsi antara dua kelompok

Estimasi proporsi dalam masing - masing kelompok diberikan oleh :

\[ \hat{p}_1 = \frac{a}{a + b} \]

\[ \hat{p}_2 = \frac{c}{c + d} \]

Estimasi proporsi gabungan :

\[ \hat{p} = \frac{x_1 + x_2}{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)}} \]

Setelah nilai Z ditemukan, kita dapat melihat keputusan apakah tolak H0 atau terima H0 dengan cara apabila |Z| lebih besar dari nilai kritis tertentu untuk tingkat signifikansi tertentu, maka H0 ditolak.

Penerapan dalam Studi Kasus

Tabel berikut merupakan data studi kasus tentang hubungan metode belajar dengan hasil tes siswa di suatu sekolah.

Lulus (+) Tidak Lulus (-) Total
Bimbingan Belajar 80 20 100
Belajar Mandiri 50 50 100
Total 130 70 200

Perhitungan Manual

Langkah 1 : Hitung Proporsi Masing-Masing Kelompok

Proporsi pada Bimbingan Belajar :

\[ \hat{p}_1 = \frac{80}{100} = 0.8 \]

Proporsi Pada Belajar Mandiri :

\[ \hat{p}_2 = \frac{50}{100} = 0.5 \]

Langkah 2 : Hitung Proporsi Gabungan

\[ \hat{p} = \frac{80 + 50}{100 + 100} = \frac{130}{200} = 0.65 \]

Langkah 3 : Hitung Statistik Uji Z

\[ Z = \frac{0.8 - 0.3}{\sqrt{0.65(1 - 0.65)\left(\frac{1}{100} + \frac{1}{100} \right)}} = \frac{0.3}{\sqrt{0.2275 \cdot 0.02}} \approx \frac{0.3}{0.0675} \approx 4.4474 \]

Perhitungan dengan R

# Data
x1 <- 80  # bimbingan belajar lulus
n1 <- 100  # total bimbingan belajar
x2 <- 50  # belajar mandiri sukses
n2 <- 100  # total belajar mandiri

# Proporsi masing-masing kelompok
p1 <- x1 / n1
p2 <- x2 / n2

# Proporsi gabungan
p_pool <- (x1 + x2) / (n1 + n2)

# Statistik uji Z
Z <- (p1 - p2) / sqrt(p_pool * (1 - p_pool) * (1/n1 + 1/n2))

# Tampilkan hasil
list(proporsi_1 = p1, proporsi_2 = p2, proporsi_gabungan = p_pool, Z = Z) 
## $proporsi_1
## [1] 0.8
## 
## $proporsi_2
## [1] 0.5
## 
## $proporsi_gabungan
## [1] 0.65
## 
## $Z
## [1] 4.447496

Interpretasi : Dengan taraf signifikansinya 5%, kita menolak H0 karena nilai |Z| yaitu 4,4474 > dari 1.96 yang berarti bahwa terdapat perbedaan proporsi antar kelompok.

Uji Asosiasi

Uji asosiasi dalam tabel kontingensi 2 × 2 bertujuan untuk mengukur hubungan antara dua variabel kategori. Untuk setiap uji asosiasi, hipotesis yang diuji adalah :

  • Hipotesis nol (H0) : tidak ada asosiasi antara dua variabel.
  • Hipotesis altenatif (H1) : terdapat asosiasi antara dua variabel.

Uji asosisasi memiliki tiga ukuran utama yaitu risk difference, relative risk, dan odds ratio.

Risk Difference (RD)

Risk difference mengukur perbedaan absolut dalam probabilitas kejadian antara dua kelompok.

Rumus : \[ RD = \frac{a}{a + b} - \frac{c}{c + d} \]Standar Error untuk RD :

\[ SE(RD) = \sqrt{ \frac{p_1(1 - p_1)}{n_1} + \frac{p_2(1 - p_2)}{n_2} } \]

Statistik Uji Z untuk RD :\[ Z = \frac{RD}{SE(RD)} \]

Relative Risk (RR)

Perbandingan antara risiko kejadian di kelompok 1 dibanding kelompok 2.

Rumus:

\[ RR = \frac{a / (a + b)}{c / (c + d)} \]

Standar Error untuk log(RR) :

\[ SE(\log(RR)) = \sqrt{ \frac{1 - p_1}{n_1 p_1} + \frac{1 - p_2}{n_2 p_2} } \]

Statistik Uji Z untuk RR :

\[ Z = \frac{\log(RR)}{SE(\log(RR))} \]

Odds Ratio (OR)

Perbandingan peluang kejadian antara dua kelompok.

Rumus: \[ OR = \frac{a \cdot d}{b \cdot c} \]Standar error untuk log(OR) :

\[ SE(\log(OR)) = \sqrt{ \frac{1}{a} + \frac{1}{b} + \frac{1}{c} + \frac{1}{d} } \]

Statistik Uji Z untuk OR :

\[ Z = \frac{\log(OR)}{SE(\log(OR))} \]

Penerapan dalam Studi Kasus

Tabel berikut merupakan data studi kasus tentang hubungan metode belajar dengan hasil tes siswa di suatu sekolah.

Lulus (+) Tidak Lulus (-) Total
Bimbingan Belajar 80 20 100
Belajar Mandiri 50 50 100
Total 130 70 200

Perhitungan Manual

  • Risk Difference :

\[ \hat{p}_1 = \frac{80}{100} = 0.8, \quad \hat{p}_2 = \frac{50}{100} = 0.5 \]

\[ RD = 0.8 - 0.5 = 0.3 \]

Standar Error RD :

\[ SE(RD) = \sqrt{ \frac{0.8 \times (1 - 0.8)}{100} + \frac{0.5 \times (1 - 0.5)}{100} }\approx 0.064 \]

Statistik Uji Z RD :

\[ Z_{RD} = \frac{0.3}{0.0.064} \approx 4.6852 \]

  • Relative Risk :

\[ RR = \frac{0.8}{0.5} = 1.6 \]

Standar Error log(RR) :

\[ SE(\log(RR)) = \sqrt{ \left( \frac{1}{80} - \frac{1}{100} \right) + \left( \frac{1}{50} - \frac{1}{100} \right) }\approx 0.1118 \]

Statistik Uji Z RR :

\[ Z_{RR} = \frac{\log(1.6)}{0.1118} \approx \frac{4.899}{0.1118} \approx 4.2038 \]

  • Odds Ratio

\[ OR = \frac{80 \cdot 50}{20 \cdot 50} = \frac{4000}{1000} = 4 \]

Standar Error log(OR) :

\[ SE(\log(OR)) = \sqrt{ \frac{1}{21} + \frac{1}{2} + \frac{1}{15} + \frac{1}{3} } \approx 0.3202 \]

Statistik Uji Z OR :

\[ Z_{OR} = \frac{\log(4)}{0.3202} \approx \frac{1.6609}{0.3202} \approx 4.33 \]

Perhitungan dengan R

# Data 
a <- 80  # Bimbingan Belajar, Lulus 
b <- 20   # Bimbingan Belajar, Tidak Lulus
c <- 50  # Belajar Mandiri, Lulus
d <- 50  # Belajar Mandiri, Tidak Lulus

# Total masing-masing grup
n1 <- a + b  # Bimbingan Belajar
n2 <- c + d  # Belajar Mandiri

# Proporsi
p1 <- a / n1
p2 <- c / n2

### --- 1. Risk Difference (RD) ---
RD <- p1 - p2

# Standard Error RD
SE_RD <- sqrt((p1 * (1 - p1)) / n1 + (p2 * (1 - p2)) / n2)

# Z-test RD
Z_RD <- RD / SE_RD

### --- 2. Relative Risk (RR) ---
RR <- p1 / p2

# SE log(RR)
SE_log_RR <- sqrt((1/a - 1/n1) + (1/c - 1/n2))

# Z-test RR
Z_RR <- log(RR) / SE_log_RR

### --- 3. Odds Ratio (OR) ---
OR <- (a * d) / (b * c)

# SE log(OR)
SE_log_OR <- sqrt(1/a + 1/b + 1/c + 1/d)

# Z-test OR
Z_OR <- log(OR) / SE_log_OR

list(Risk_Difference = RD, SE_RD = SE_RD, Z_RD = Z_RD, Relative_Risk = RR, SE_log_RR = SE_log_RR, Z_RR = Z_RR, Odds_Ratio = OR, SE_log_OR = SE_log_OR, Z_OR = Z_OR)
## $Risk_Difference
## [1] 0.3
## 
## $SE_RD
## [1] 0.06403124
## 
## $Z_RD
## [1] 4.685213
## 
## $Relative_Risk
## [1] 1.6
## 
## $SE_log_RR
## [1] 0.1118034
## 
## $Z_RR
## [1] 4.20384
## 
## $Odds_Ratio
## [1] 4
## 
## $SE_log_OR
## [1] 0.3201562
## 
## $Z_OR
## [1] 4.330056

Uji Independensi

Uji independensi digunakan untuk menentukan apakah ada hubungan statistik antara dua variabel kategorikal.

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} \]

  • \(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.
Contoh perhitungan

Tabel berikut merupakan data studi kasus tentang hubungan metode belajar dengan hasil tes siswa di suatu sekolah.

Lulus (+) Tidak Lulus (-) Total
Bimbingan Belajar 80 20 100
Belajar Mandiri 50 50 100
Total 130 70 200

Perhitungan Manual :

Langkah 1 : Hitung nilai yang diharapkan (E)

Nilai yang diharapkan \(E_{ij}\) dihitung dengan rumus :

\[ E_{ij} = \frac{R_i \times C_j}{N} \]

  • \(R_i\) adalah total baris ke-\(i\)
  • \(C_j\) adalah total kolom ke-\(j\)
  • \(N\) adalah total sampel

Menghitung nilai yang diharapkan:

Untuk Bimbingan Belajar, Lulus (\(E_{11}\)):

\[ E_{11} = \frac{(100 \times 130)}{200} = 65 \]

Untuk Bimbingan Belajar, Tidak Lulus (\(E_{12}\)):

\[ E_{12} = \frac{(100 \times 70)}{200} = 35 \]

Untuk Belajar Mandiri, Lulus (\(E_{21}\)):

\[ E_{21} = \frac{(100 \times 130)}{200} = 65 \]

Untuk Belajar Mandiri, Tidak Lulus (\(E_{22}\)):

\[ E_{22} = \frac{(100 \times 70)}{200} = 35 \]

Langkah 2 : Hitung Chi-Square (\(\chi^2\))

\[ \chi^2 = \sum \frac{(O - E)^2}{E} \]

  • \(O\) adalah nilai observasi
  • \(E\) adalah nilai yang diharapkan

Untuk Bimbingan Belajar, Lulus (\(O_{11} = 80\), \(E_{11} = 65\)):

\[ \frac{(80 - 65)^2}{65} = \frac{225}{65} = 3.4615 \]

Untuk Bimbingan Belajar, Tidak Lulus (\(O_{12} = 20\), \(E_{12} = 35\)):

\[ \frac{(20 - 35)^2}{35} = \frac{225}{35} = 6.4286 \]

Untuk Belajar Mandiri, Lulus (\(O_{21} = 50\), \(E_{21} = 65\)):

\[ \frac{(50 - 65)^2}{65} = \frac{225}{65} = 3.4615 \]

Untuk Belajar Mandiri, Tidak Lulus (\(O_{22} = 50\), \(E_{22} = 35\)):

\[ \frac{(50 - 35)^2}{35} = \frac{225}{35} = 6.4286 \]

Chi-Square Total :

\[ \chi^2 = 3.4615 + 6.4286 + 3.4615 + 6.4286 = 18.484 \]

Perhitungan dengan R :

# Membuat data tabel 2x2
data <- matrix(c(80, 20, 50, 50), nrow = 2, byrow = TRUE)

# Menambahkan nama baris dan kolom
rownames(data) <- c("Bimbingan Belajar", "Belajar Mandiri")
colnames(data) <- c("Lulus", "Tidak Lulus")

# Menampilkan data tabel
data
##                   Lulus Tidak Lulus
## Bimbingan Belajar    80          20
## Belajar Mandiri      50          50
# Menghitung chi-square
chi_square_test <- chisq.test(data)

# Menampilkan hasil uji chi-square
chi_square_test
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data
## X-squared = 18.484, df = 1, p-value = 1.714e-05
list(data = data, chi_square_test = chi_square_test)
## $data
##                   Lulus Tidak Lulus
## Bimbingan Belajar    80          20
## Belajar Mandiri      50          50
## 
## $chi_square_test
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data
## X-squared = 18.484, df = 1, p-value = 1.714e-05

Interpretasi : Dengan taraf signifkansi 5% didapatkan p-value < alpha, maka tolak H0 yang berarti metode belajar (bimbingan belajar vs belajar mandiri) berpengaruh terhadap kemungkinan lulus atau tidak lulus siswa.

Uji Partisi Chi Square

Partisi Chi-Square Partisi Chi-Square digunakan untuk mengidentifikasi kategori mana dalam tabel kontingensi yang bertanggung jawab atas hubungan yang signifikan. Jika uji Chi-Square pada tabel kontingensi I×J signifikan, maka partisi Chi-Square memungkinkan kita untuk menguraikan efek hubungan dalam subkelompok yang lebih kecil.

Contoh Perhitungan

Tabel berikut merupakan data studi kasus tentang hubungan konsumsi minuman soda (botol per hari) dengan risiko terkena diabetes.

0 botol / hari 1-3 botol / hari >3 botol / hari Total
Kontrol 25 25 12 62
Diabetes 0 1 3 4
Total 25 26 15 66

Perhitungan Manual

Langkah 1 : Menghitung Frekuensi yang Diharapkan

Frekuensi yang diharapkan dihitung menggunakan rumus:

\[ E_i = \frac{(Total \, \text{Baris}) \times (Total \, \text{Kolom})}{Grand \, \text{Total}} \]

Berdasarkan tabel di atas, kita dapat menghitung frekuensi yang diharapkan untuk setiap sel sebagai berikut:

\[ E_{1,1} = \frac{62 \times 25}{66} = 23.636 \]

\[ E_{1,2} = \frac{62 \times 26}{66} = 24.545 \]

\[ E_{1,3} = \frac{62 \times 15}{66} = 14.091 \]

\[ E_{2,1} = \frac{4 \times 25}{66} = 1.515 \]

\[ E_{2,2} = \frac{4 \times 26}{66} = 1.515 \]

\[ E_{2,3} = \frac{4 \times 15}{66} = 0.909 \]

Langkah 2: Menghitung Statistik Chi-Square

Sekarang kita dapat menghitung statistik chi-square (\(\chi^2\)) dengan rumus:

\[ \chi^2 = \sum \frac{(O_i - E_i)^2}{E_i} \]

Dimana \(O_i\) adalah frekuensi yang diamati dan \(E_i\) adalah frekuensi yang diharapkan. Berdasarkan data yang ada, kita menghitung chi-square untuk setiap sel:

\[ \chi^2 = \frac{(25 - 23.636)^2}{23.636} + \frac{(25 - 24.545)^2}{24.545} + \frac{(12 - 14.091)^2}{14.091} + \frac{(0 - 1.515)^2}{1.515} + \frac{(1 - 1.515)^2}{1.515} + \frac{(3 - 0.909)^2}{0.909} \]

\[ \chi^2 = 0.0785 + 0.0084 + 0.3107 + 1.513 + 0.1754 + 4.809 = 6.9562 \]Pehitungan Dengan R :

# Membuat matriks data
data_matrix <- matrix(c(25, 25, 12, 0, 1, 3), nrow = 2, byrow = TRUE)
colnames(data_matrix) <- c("0 botol / hari", "1-3 botol / hari", ">3 botol / hari")
rownames(data_matrix) <- c("Kontrol", "Diabetes")

# Melakukan uji chi-square
chi_square_test <- chisq.test(data_matrix)
## Warning in chisq.test(data_matrix): Chi-squared approximation may be incorrect
# Menampilkan hasil uji
chi_square_test
## 
##  Pearson's Chi-squared test
## 
## data:  data_matrix
## X-squared = 6.9562, df = 2, p-value = 0.03087

Intepretasi : Dengan taraf signifikansi 5%, maka kita menolak H0 dikarenakan p value (0.03087) lebih kecil dari alpha (0.05). Dikarenakan menolak H0, maka terdapat hubungan signifikan antara konsumsi minuman soda dan risiko diabetes.

Uji Likelihood Ratio

Uji Likelihood Ratio (G²) adalah alternatif dari uji chi-square yang digunakan untuk menguji hipotesis independensi dalam tabel kontingensi \(I \times J\). Statistik uji ini diberikan oleh :

\[ G^2 = 2 \sum_{i} \sum_{j} n_{ij} \ln \left( \frac{n_{ij}}{\hat{\mu}_{ij}} \right) \]

  • \(n_{ij}\) adalah frekuensi observasi dalam tabel kontingensi.

  • \(\hat{\mu}_{ij} = n \cdot p_i \cdot p_j\) adalah frekuensi yang diharapkan.

  • \(p_i\) adalah proporsi untuk baris \(i\).

  • \(p_j\) adalah proporsi untuk kolom \(j\).

  • \(n\) adalah total frekuensi.
Contoh Perhitungan

Tabel berikut merupakan data studi kasus tentang hubungan metode belajar dengan hasil tes siswa di suatu sekolah.

Lulus (+) Tidak Lulus (-) Total
Bimbingan Belajar 80 20 100
Belajar Mandiri 50 50 100
Total 130 70 200

Langkah 1: Hitung Frekuensi Ekspektasi

Gunakan rumus: \[ \hat{\mu}_{ij} = \frac{(\text{Total Baris}_i) \times (\text{Total Kolom}_j)}{\text{Grand Total}} \]

\[ \begin{aligned} \hat{\mu}_{11} &= \frac{100 \times 130}{200} = 65 \\ \hat{\mu}_{12} &= \frac{100 \times 70}{200} = 35 \\ \hat{\mu}_{21} &= \frac{100 \times 130}{200} = 65 \\ \hat{\mu}_{22} &= \frac{100 \times 70}{200} = 35 \\ \end{aligned} \]

Langkah 2: Hitung Statistik Uji G²

\[ \begin{aligned} G^2 &= 2 \sum O_{ij} \ln \left( \frac{O_{ij}}{E_{ij}} \right) \\ &= 2 \left[ 80 \ln \left( \frac{80}{65} \right) + 20 \ln \left( \frac{20}{35} \right) + 50 \ln \left( \frac{50}{65} \right) + 50 \ln \left( \frac{50}{35} \right) \right] \\ &= 20.269 \end{aligned} \]

Langkah 3 : Bandingkan dengan Distribusi Chi-Square

Derajat bebas: \[ (I - 1)(J - 1) = (2 - 1)(2 - 1) = 1 \]

Nilai kritis \(\chi^2\) untuk \(\alpha = 0.05\) dan df = 1 adalah 3.841.

Karena \(G^2 = 20.269 > 3.841\), maka tolak H0, artinya terdapat cukup bukti yang menunjukkan bahwa metode belajar berpengaruh signifikan terhadap kelulusan siswa.

Perhitungan dengan R :

# Data Observasi
observed <- matrix(c(80, 20,
                     50, 50), 
                   nrow = 2, byrow = TRUE)

# Menambahkan nama baris dan kolom
rownames(observed) <- c("Bimbingan Belajar", "Belajar Mandiri")
colnames(observed) <- c("Lulus", "Tidak Lulus")

# Total per baris dan kolom
row_totals <- rowSums(observed)
col_totals <- colSums(observed)
grand_total <- sum(observed)

# Frekuensi harapan (Expected frequencies)
expected <- outer(row_totals, col_totals) / grand_total

# Hitung G²
G2 <- 2 * sum(observed * log(observed / expected), na.rm = TRUE)

# Derajat bebas
df <- (nrow(observed) - 1) * (ncol(observed) - 1)

# Nilai kritis chi-square
chi_critical <- qchisq(0.95, df)

# Output hasil
list(statistik_G = G2, Derajat_bebas = df, Nilai_kritis = chi_critical)
## $statistik_G
## [1] 20.26873
## 
## $Derajat_bebas
## [1] 1
## 
## $Nilai_kritis
## [1] 3.841459
if (G2 > chi_critical) {
  cat("Kesimpulan: Tolak H0 - Ada hubungan signifikan.\n")
} else {
  cat("Kesimpulan: Gagal tolak H0 - Tidak ada hubungan signifikan.\n")
}
## Kesimpulan: Tolak H0 - Ada hubungan signifikan.

Tabel Kontingensi 3 Arah

Tabel kontingensi tiga arah merupakan perluasan dari tabel kontingensi dua arah yang digunakan untuk menganalisis hubungan antara tiga variabel kategori secara bersamaan. Dalam berbagai situasi, keterkaitan antara dua variabel (misalnya X dan Y) dapat dipengaruhi oleh variabel ketiga Z, yang dikenal sebagai variabel kontrol atau kovariat. Tabel ini sering dimanfaatkan dalam analisis data ketika terdapat kemungkinan adanya faktor pengganggu yang dapat memengaruhi hubungan antara dua variabel utama. Contohnya :

Fungsi Tabel Marginal

Tabel Parsial dan Marginal

Tabel parsial merupakan tabel yang menyajikan hubungan antara X dan Y berdasarkan tiap kategori Z, sedangkan tabel marginal adalah tabel yang mengabaikan Z, dengan menjumlahkan data dari seluruh kategori Z.

Contoh :

Tabel berikut menunjukkan data mengenai peristiwa cedera pada kecelakaan lalu lintas berdasarkan jenis kendaraan dan waktu berkendara.

Waktu Berkendara Jenis Kendaraan Cedera Tidak Cedera
Malam Motor 43 134
Malam Mobil 26 149
Siang Motor 29 23
Siang Mobil 22 26

Tabel Parsial

Tabel frekuensi parsial menyajikan hubungan antara dua variabel kategori dalam tabel kontingensi tiga arah dengan mempertimbangkan satu variabel sebagai kontrol. Tabel ini membantu dalam menambah hubungan bersyarat antara variabel dalam analisis data kategori.

  • Tabel Frekuensi Parsial untuk Z (Waktu Berkendara) = Malam : \[ \begin{array}{|c|c|c|} \hline \text{Jenis Kendaraan} & \text{Cedera} & \text{Tidak Cedera} \\ \hline \text{Motor} & 43 & 134 \\ \text{Mobil} & 26 & 149 \\ \hline \end{array} \]

Tabel Frekuensi Parsial untuk Z (Waktu Berkendara) = Malam menunjukkan hubungan antara jenis kendaraan dan apakah individu tersebut cedera atau tidak saat kecelakaan di malam hari.

  • Tabel Frekuensi Parsial untuk Z (Waktu Berkendara) = Siang : \[ \begin{array}{|c|c|c|} \hline \text{Jenis Kendaraan} & \text{Cedera} & \text{Tidak Cedera} \\ \hline \text{Motor} & 29 & 23 \\ \text{Mobil} & 22 & 26 \\ \hline \end{array} \]

Tabel Frekuensi Parsial untuk Z (Waktu Berkendara) = Siang menunjukkan hubungan antara jenis kendaraan dan apakah individu tersebut cedera atau tidak saat kecelakaan di siang hari.

Dengan R :

# Membuat data dengan array
data3 <- array(
  c(43, 26, 134, 149, 29, 22, 23, 36), 
  dim = c(2, 2, 2), 
  dimnames = list(
    Jenis_Kendaraan = c("Motor", "Mobil"),
    Intercourse = c("Cedera", "Tidak Cedera"),
    Ras = c("Malam", "Siang")
  )
)

# Tampilkan array untuk memastikan datanya
data3
## , , Ras = Malam
## 
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     43          134
##           Mobil     26          149
## 
## , , Ras = Siang
## 
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     29           23
##           Mobil     22           36
# Tabel parsial: hanya untuk waktu malam
parsial_malam <- data3[, , "Malam"]
parsial_malam
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     43          134
##           Mobil     26          149
# Tabel parsial: hanya untuk waktu siang
parsial_siang <- data3[, , "Siang"]
parsial_siang
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     29           23
##           Mobil     22           36
list(parsial_malam = parsial_malam, parsial_siang = parsial_siang)
## $parsial_malam
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     43          134
##           Mobil     26          149
## 
## $parsial_siang
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     29           23
##           Mobil     22           36

Tabel Marginal

Tabel frekuensi marginal menampilkan jumlah total observasi untuk setiap variabel dengan mengabaikan variabel lainnya dalam tabel kontingensi tiga arah. Tabel ini membantu dalam memahami distribusi kategori secara agregat tanpa mempertimbangkan hubungan antarvariabel. Tabel frekuensi marginal dihitung dengan menjumlahkan frekuensi dari tabel kontingensi tiga arah berdasarkan variabel yang tersisa.

  • Tabel Frekuensi Marginal untuk waktu berkendara dan peristiwa cedera
Waktu Berkendara Cedera Tidak Cedera Total
Malam 69 283 352
Siang 51 49 100
  • Tabel Frekuensi Marginal untuk jenis kendaraan dan peristiwa cedera
Jenis Kendaraan Cedera Tidak Cedera Total
Motor 72 229 301
Mobil 48 175 223

Dengan R :

# Membuat data dengan array
data3 <- array(
  c(43, 26, 134, 149, 29, 22, 23, 36), 
  dim = c(2, 2, 2), 
  dimnames = list(
    Jenis_Kendaraan = c("Motor", "Mobil"),
    Intercourse = c("Cedera", "Tidak Cedera"),
    Ras = c("Malam", "Siang")
  )
)

# Tampilkan array untuk memastikan datanya
data3
## , , Ras = Malam
## 
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     43          134
##           Mobil     26          149
## 
## , , Ras = Siang
## 
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     29           23
##           Mobil     22           36
# Menjumlahkan berdasarkan waktu
tabel_marginal_waktu <- apply(data3, c(1, 2), sum)
tabel_marginal_waktu
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     72          157
##           Mobil     48          185
# Menjumlahkan berdasarkan jenis kendaraan
tabel_marginal_kendaraan <- apply(data3, c(3, 2), sum)
tabel_marginal_kendaraan
##        Intercourse
## Ras     Cedera Tidak Cedera
##   Malam     69          283
##   Siang     51           59
list(tabel_marginal_waktu = tabel_marginal_waktu, tabel_marginal_kendaraan = tabel_marginal_kendaraan)
## $tabel_marginal_waktu
##                Intercourse
## Jenis_Kendaraan Cedera Tidak Cedera
##           Motor     72          157
##           Mobil     48          185
## 
## $tabel_marginal_kendaraan
##        Intercourse
## Ras     Cedera Tidak Cedera
##   Malam     69          283
##   Siang     51           59

Distribusi Peluang

Distribusi peluang pada tabel kontingensi tiga arah adalah representasi probabilitas relatif dari setiap kombinasi tiga variabel kategorikal dalam tabel. Alih-alih menunjukkan jumlah kasus atau frekuensi, distribusi peluang menunjukkan seberapa besar kemungkinan suatu kombinasi kategori terjadi dibandingkan dengan seluruh data. Distribusi ini sangat berguna untuk menganalisis hubungan antar variabel dalam skala probabilistik, baik secara keseluruhan, sebagian, maupun dalam konteks tertentu. Distribusi pada tabel kontingensi tiga arah terdira dari peluang bersama, peluang marginal, dan peluang bersyarat.

Peluang Bersama

Peluang bersama didefinisikan sebagai:

\[ P(Z, X, Y) = \frac{f(Z, X, Y)}{N} \]

Perhitungan dengan R :

# Data awal
data <- data.frame(
  Waktu = c("Malam", "Malam", "Siang", "Siang"),
  Kendaraan = c("Motor", "Mobil", "Motor", "Mobil"),
  Cedera = c(43, 26, 29, 22),
  Tidak_Cedera = c(134, 149, 23, 26)
)

# Total seluruh kejadian
total <- sum(data$Cedera + data$Tidak_Cedera)

# Peluang bersama untuk masing-masing kondisi
data$Peluang_Cedera <- data$Cedera / total
data$Peluang_Tidak_Cedera <- data$Tidak_Cedera / total

# Tampilkan peluang bersama
cat("=== Peluang Bersama ===\n")
## === Peluang Bersama ===
print(data[, c("Waktu", "Kendaraan", "Peluang_Cedera", "Peluang_Tidak_Cedera")])
##   Waktu Kendaraan Peluang_Cedera Peluang_Tidak_Cedera
## 1 Malam     Motor     0.09513274           0.29646018
## 2 Malam     Mobil     0.05752212           0.32964602
## 3 Siang     Motor     0.06415929           0.05088496
## 4 Siang     Mobil     0.04867257           0.05752212

Interpretasi :

  • Peluang bahwa pengendara motor malam hari mengalami cedera adalah 11.00%.

  • Peluang bahwa pengendara mobil siang hari tidak mengalami cedera adalah 6.65%.

Peluang Marginal

Peluang marginal didefinisikan sebagai berikut

  • untuk waktu berkendara : \[ P(Z) = \frac{n(Z)}{N} \]

Perhitungan dengan R :

# Peluang marginal berdasarkan waktu berkendara
marginal_waktu <- aggregate(cbind(Cedera, Tidak_Cedera) ~ Waktu, data = data, sum)
marginal_waktu$Peluang_Marginal <- (marginal_waktu$Cedera + marginal_waktu$Tidak_Cedera) / total

# Tampilkan hasil
print(marginal_waktu[, c("Waktu", "Peluang_Marginal")])
##   Waktu Peluang_Marginal
## 1 Malam        0.7787611
## 2 Siang        0.2212389

Interpretasi : 90.06% dari semua kejadian dalam data terjadi saat malam hari, artinya mayoritas insiden terjadi pada malam hari.

  • untuk jenis kendaraan : \[ P(X) = \frac{n(X)}{N} \]

Perhitungan dengan R :

# Peluang marginal berdasarkan jenis kendaraan
marginal_kendaraan <- aggregate(cbind(Cedera, Tidak_Cedera) ~ Kendaraan, data = data, sum)
marginal_kendaraan$Peluang_Marginal <- (marginal_kendaraan$Cedera + marginal_kendaraan$Tidak_Cedera) / total

cat("\n=== Peluang Marginal Berdasarkan Kendaraan ===\n")
## 
## === Peluang Marginal Berdasarkan Kendaraan ===
print(marginal_kendaraan[, c("Kendaraan", "Peluang_Marginal")])
##   Kendaraan Peluang_Marginal
## 1     Mobil        0.4933628
## 2     Motor        0.5066372

Interpretasi : 58.61% dari total kejadian melibatkan kendaraan motor, menunjukkan bahwa motor lebih banyak terlibat dibanding mobil.

Peluang Bersyarat

Mencari peluang seseorang berkendara pada malam hari dan mengalami cedera : \[ P(Cedera| Waktu = Malam) = \frac{P(Cedera,Malam)}{P(Malam)} \]

# P(Cedera | Waktu = Malam)
malam <- subset(data, Waktu == "Malam")
P_cedera_given_malam <- sum(malam$Cedera) / sum(malam$Cedera + malam$Tidak_Cedera)

#Tampilkan hasil
cat("=== Peluang Bersyarat ===\n")
## === Peluang Bersyarat ===
cat("P(Cedera | Waktu = Malam):", round(P_cedera_given_malam, 4), "\n")
## P(Cedera | Waktu = Malam): 0.196

Interpretasi : Dapat diketahui bahwa berkendara di malam hari, maka peluang mereka mengalami cedera adalah 13.73%.

Ukuran Asosiasi

Perbedaan Peluang (Risk Difference, RD)

Risk Difference didefinisikan sebagai :

\[ RD = \frac{a}{a + b} - \frac{c}{c + d} \]

Perhitungan dengan R :

# Risk Difference
cedera_malam <- 69
total_malam <- 352
risk_malam <- cedera_malam / total_malam

cedera_siang <- 51
total_siang <- 100
risk_siang <- cedera_siang / total_siang

risk_difference <- risk_malam - risk_siang
cat("Risk Difference (Malam - Siang):", round(risk_difference, 4), "\n")
## Risk Difference (Malam - Siang): -0.314

Interpretasi :

  • RD negatif menunjukkan siang hari lebih berisiko.

  • Jika hasil = -0.314, maka risiko cedera berkurang 31.4% saat malam dibanding siang.

Risiko Relatif (Relative Risk, RR)

Relative Risk didefinisikan sebagai :

\[ RR = \frac{a / (a + b)}{c / (c + d)} \]

Perhitungan dengan R :

# Relative Risk
relative_risk <- risk_malam / risk_siang
cat("Relative Risk (Malam vs Siang):", round(relative_risk, 4), "\n")
## Relative Risk (Malam vs Siang): 0.3844

Interpretasi :

  • RR < 1 menunjukkan risiko lebih tinggi di siang hari.

  • Jika RR = 0.3843, maka risiko cedera malam hari hanya 38.4% dari siang hari.

Odds Ratio (OR)

Odds Ratio didefinisikan sebagai :

\[ OR = \frac{a \cdot d}{b \cdot c} \]

Perhitungan dengan R :

# Odds Ratio
tidak_cedera_malam <- 283
tidak_cedera_siang <- 49

odds_malam <- cedera_malam / tidak_cedera_malam
odds_siang <- cedera_siang / tidak_cedera_siang

odds_ratio <- odds_malam / odds_siang
cat("Odds Ratio (Malam vs Siang):", round(odds_ratio, 4), "\n")
## Odds Ratio (Malam vs Siang): 0.2343

Interpretasi :

  • OR < 1 menunjukkan peluang cedera lebih kecil di malam hari.

  • Jika OR = 0.2378, maka peluang cedera sekitar 76% lebih rendah saat malam dibanding siang.

Inferensi Tabel Kontingensi Tiga Arah

Tabel kontingensi tiga arah digunakan untuk mengevaluasi hubungan antara dua variabel kategorik, misalnya adalah hubungan antara jenis kendaraan (X) dan risiko cedera saat kecelakaan (Y) dengan variabel waktu berkendara (Z).. Dalam analisis ini, dibuat beberapa tabel 2×2 untuk setiap kategori variabel ketiga serta satu tabel gabungan (marginal) yang mengabaikan variabel tersebut.

Jika nilai odds ratio dari setiap tabel parsial relatif stabil, maka kita dapat menghitung odds ratio gabungan menggunakan metode Mantel-Haenszel.

Independensi Bersyarat dalam Tabel Kontingensi Tiga Arah

Independensi bersyarat berarti bahwa dua variabel (misalnya X dan Y) tidak saling berhubungan pada setiap level dari variabel ketiga (Z). Dengan kata lain, X dan Y dikatakan independen secara bersyarat terhadap Z jika odds ratio di masing-masing strata Z bernilai 1.

Pengujian Statistik untuk Independensi Bersyarat

Untuk menguji apakah X dan Y benar-benar independen bersyarat terhadap Z, digunakan uji statistik Cochran-Mantel-Haenszel (CMH). Uji ini berguna untuk melihat hubungan antara dua variabel kategori sambil memperhitungkan potensi variabel perancu (confounder).

\[ CMH = \frac{\sum_k \left( n_{1ik} - \mu_{1ik} \right)^2}{\sum_k \text{var}(n_{1ik})} \]

  • \(n_{1ik}\): nilai frekuensi sel baris 1 kolom 1 pada tabel parsial ke-\(k\).
  • \(\mu_{1ik}\): nilai ekspektasi sel baris 1 kolom 1 pada tabel parsial ke-\(k\), dihitung dengan rumus:

\[ \mu_{1ik} = E(n_{1ik}) = \frac{n_{1.k} \cdot n_{1ik}}{n_{.k}} \]

Varians dari \(n_{1ik}\) diberikan oleh :

\[ \text{var}(n_{1ik}) = \frac{n_{1.k} \cdot n_{2.k} \cdot n_{1ik} \cdot n_{2.k}}{n_{2.k}^2 \cdot (n_{.k} - 1)} \]

Penerapan dalam Studi Kasus

Tabel berikut menunjukkan data mengenai peristiwa cedera pada kecelakaan lalu lintas berdasarkan jenis kendaraan dan waktu berkendara.

Waktu Berkendara Jenis Kendaraan Cedera Tidak Cedera
Malam Motor 43 134
Malam Mobil 26 149
Siang Motor 29 23
Siang Mobil 22 26

Perhitungan Manual :

  1. Rumus Nilai Harapan \(\mu_{1ik}\) untuk Setiap Strata

Rumus nilai harapan \(\mu_{1ik}\) untuk setiap strata \(k\):

\[ \mu_{1ik} = \frac{(n_{i+k} \times n_{+ik})}{n_{++k}} \]

  • \(n_{i+k}\) = Total laki - laki dalam strata \(k\) - \(n_{+ik}\) = Total Intercourse dalam strata \(k\) - \(n_{++k}\) = Total individu dalam strata \(k\)

Untuk Malam (k=1) :

\[ \mu_{111} = \frac{(177 \times 69)}{352} = 34.696 \]

Untuk Siang (k=2) :

\[ \mu_{112} = \frac{(52 \times 51)}{110} = 24,109 \]

  1. Menghitung Varians \(\text{Var}(n_{1ik})\)

Rumus varians:

\[ \text{Var}(n_{1ik}) = \frac{n_{i+k} \times n_{0+k} \times n_{+ik} \times n_{+0k}}{n_{++k}^2 \times (n_{++k} - 1)} \]

Untuk Malam (k=1) :

\[ \text{Var}(n_{111}) = \frac{(177 \times 134 \times 43 \times 149)}{352^2 \times (352 - 1)} = 13.91 \]

Untuk Siang (k=2) :

\[ \text{Var}(n_{112}) = \frac{(175 \times 149 \times 26 \times 134)}{350^2 \times (350 - 1)} = 9.91 \]

  1. Menghitung Statistik CMH

Rumus statistik CMH:

\[ X^2_{CMH} = \frac{\sum_k (n_{1ik} - \mu_{1ik})^2}{\sum_k \text{Var}(n_{1ik})} \]

\[ X^2_{CMH} = \frac{((43 - 34.696) + (29 - 24.109))^2}{13.91 + 9.91} = 7.309 \]

Intepretasi : Dengan taraf signifikansi 5%, maka kita memiliki nilai kritis sebesar 3.841. Dikarenakan 7.309 > 3.841 maka hipotesis nol ditolak yang berarti terdapat hubungan signifikan antara jenis kendaraan terhadap risiko cedera setelah mengontrol waktu berkendara.

Perhitungan dengan R :

#Membuat data
data_cmh <- array(
  c(43, 26, 134, 149,  # Malam
    29, 22, 23, 36), # Siang
  dim = c(2, 2, 2),
  dimnames = list(
    "Jenis Kendaraan" = c("Motor", "Mobil"),
    "Risiko Cedera" = c("Cedera", "Tidak Cedera"),
    "Waktu Berkendara" = c("Malam", "Siang")
  )
)

# Menampilkan tabel
print(data_cmh)
## , , Waktu Berkendara = Malam
## 
##                Risiko Cedera
## Jenis Kendaraan Cedera Tidak Cedera
##           Motor     43          134
##           Mobil     26          149
## 
## , , Waktu Berkendara = Siang
## 
##                Risiko Cedera
## Jenis Kendaraan Cedera Tidak Cedera
##           Motor     29           23
##           Mobil     22           36
cmh_base <- mantelhaen.test(data_cmh, correct = FALSE)
cmh_base
## 
##  Mantel-Haenszel chi-squared test without continuity correction
## 
## data:  data_cmh
## Mantel-Haenszel X-squared = 8.3751, df = 1, p-value = 0.003804
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
##  1.229343 2.967939
## sample estimates:
## common odds ratio 
##          1.910135

Intepretasi : Dikarenakan p value 0.003804 < 0.05 maka keputusannya adalah menolak H0 yang berarti terdapat hubungan signifikan antara jenis kelamin terhadap perlakuan intercourse setelah mengontrol variabel ras.

Generalized Linear Model dalam Analisis Data Kategorik

Generalized Linear Model (GLM) adalah kerangka statistik yang memperluas regresi linier untuk menangani berbagai jenis distribusi respons, termasuk data kategorik seperti biner atau multikategori. GLM terdiri dari tiga komponen utama: distribusi respons dari keluarga eksponensial, fungsi link yang menghubungkan mean respons dengan prediktor linier, dan prediktor linier itu sendiri. Bagian ini menjelaskan secara mendetail keluarga eksponensial (dengan distribusi-distribusinya), model regresi logistik, dan model regresi Poisson, yang semuanya relevan untuk analisis data kategorik.

Exponential Family

Keluarga eksponensial adalah kelompok distribusi probabilitas yang memiliki bentuk fungsi densitas atau massa probabilitas yang seragam, memungkinkan fleksibilitas dalam pemodelan data kategorik dan lainnya. Distribusi dalam keluarga ini mencakup Bernoulli, Binomial, Multinomial, Poisson, dan Normal, yang masing-masing memiliki parameter spesifik untuk memodelkan mean dan varians respons. Bentuk umum fungsi densitas (atau massa) untuk keluarga eksponensial adalah:

\[ f(y; \theta, \phi) = \exp\left[\frac{y\theta - b(\theta)}{a(\phi)} + c(y, \phi)\right] \]

di mana:

\(y\): Variabel respons. \(\theta\): Parameter alami (natural parameter). \(\phi\): Parameter dispersi (dispersion parameter). \(b(\theta)\): Fungsi yang menentukan mean dan varians. \(a(\phi)\): Fungsi skala, sering konstan atau bergantung pada \(\phi\). \(c(y, \phi)\): Fungsi normalisasi.

Mean dan varians respons diberikan oleh:

Mean: \(E(Y) = \mu = b'(\theta)\). Varians: \(\text{Var}(Y) = b''(\theta) a(\phi)\).

Berikut adalah distribusi-distribusi utama dalam keluarga eksponensial yang relevan untuk data kategorik:

Distribusi Bernoulli

Distribusi Bernoulli memodelkan data biner dengan dua hasil (sukses = 1, gagal = 0) dalam satu percobaan. Ini adalah kasus dasar untuk data kategorik biner, seperti respons ya/tidak.

Fungsi massa probabilitasnya adalah: \[ P(Y = y) = \pi^y (1-\pi)^{1-y}, \quad y = 0, 1 \]

Dalam bentuk eksponensial: \[ f(y; \pi) = \exp\left[y \log\left(\frac{\pi}{1-\pi}\right) + \log(1-\pi)\right] \]

di mana:

\(\theta = \log\left(\frac{\pi}{1-\pi}\right)\) (logit dari probabilitas sukses). \(b(\theta) = \log(1 + e^\theta)\). \(a(\phi) = 1\) (tidak ada parameter dispersi). \(c(y, \phi) = 0\).

Parameter, Mean, dan Varians:

  • Parameter: \(\pi\) (probabilitas sukses, \(0 \leq \pi \leq 1\)).

  • Mean: \(E(Y) = \pi\).

  • Varians: \(\text{Var}(Y) = \pi(1-\pi)\).

Distribusi Bernoulli digunakan dalam regresi logistik untuk memodelkan probabilitas kejadian biner, seperti apakah pasien sembuh atau tidak setelah pengobatan.

Distribusi Binomial

Distribusi Binomial memodelkan jumlah keberhasilan dalam \(n\) percobaan independen, masing-masing dengan probabilitas sukses \(\pi\). Ini cocok untuk data kategorik biner yang diulang.

Fungsi massa probabilitasnya adalah:

\[ P(Y = y) = \binom{n}{y} \pi^y (1-\pi)^{n-y}, \quad y = 0, 1, \ldots, n \] Dalam bentuk eksponensial untuk \(Y/n\) (proporsi sukses):

\[ f(y; \pi, n) = \exp\left[y \log\left(\frac{\pi}{1-\pi}\right) + n \log(1-\pi) + \log\binom{n}{y}\right] \]

di mana:

\(\theta = \log\left(\frac{\pi}{1-\pi}\right)\). \(b(\theta) = n \log(1 + e^\theta)\). \(a(\phi) = 1\). \(c(y, \phi) = \log\binom{n}{y}\).

Parameter, Mean, dan Varians:

  • Parameter: \(n\) (jumlah percobaan), \(\pi\) (probabilitas sukses).

  • Mean: \(E(Y) = n\pi\).

  • Varians: \(\text{Var}(Y) = n\pi(1-\pi)\).

Distribusi Binomial digunakan untuk menganalisis frekuensi keberhasilan dalam percobaan berulang, seperti jumlah pasien yang sembuh dari \(n\) pasien dalam uji klinis.

Distribusi Multinomial

Distribusi Multinomial memodelkan pembagian \(n\) observasi independen ke dalam \(k\) kategori, masing-masing dengan probabilitas \(\pi_i\). Ini cocok untuk data kategorik multikategori, seperti preferensi konsumen.

Fungsi massa probabilitasnya adalah: \[ P(n_1, n_2, \ldots, n_k) = \frac{n!}{n_1! n_2! \cdots n_k!} \pi_1^{n_1} \pi_2^{n_2} \cdots \pi_k^{n_k} \]

di mana \(\sum n_i = n\) dan \(\sum \pi_i = 1\). Dalam bentuk eksponensial: \[ f(n_1, \ldots, n_k; \pi_1, \ldots, \pi_k) = \exp\left[\sum_{i=1}^{k-1} n_i \log\left(\frac{\pi_i}{\pi_k}\right) + \log\left(\frac{n!}{n_1! \cdots n_k!}\right)\right] \]

di mana:

\(\theta_i = \log\left(\frac{\pi_i}{\pi_k}\right)\) untuk kategori \(i = 1, \ldots, k-1\). \(b(\theta) = n \log\left(1 + \sum_{i=1}^{k-1} e^{\theta_i}\right)\). \(a(\phi) = 1\). \(c(n_1, \ldots, n_k) = \log\left(\frac{n!}{n_1! \cdots n_k!}\right)\).

Parameter, Mean, dan Varians:

  • Parameter: \(n\) (jumlah observasi), \(\pi_1, \ldots, \pi_k\) (probabilitas kategori).

  • Mean: \(E(n_i) = n \pi_i\).

  • Varians: \(\text{Var}(n_i) = n \pi_i (1 - \pi_i)\).

  • Kovarians: \(\text{Cov}(n_i, n_j) = -n \pi_i \pi_j\) untuk \(i \neq j\).

Distribusi Multinomial digunakan dalam regresi logistik multinomial untuk memodelkan pilihan konsumen atau klasifikasi ke dalam beberapa kategori, seperti tingkat kepuasan.

Distribusi Poisson

Distribusi Poisson memodelkan jumlah kejadian acak dalam interval waktu atau ruang tertentu, cocok untuk data kategorik berupa frekuensi atau hitungan dalam kategori.

Fungsi massa probabilitasnya adalah: \[P(Y = y) = \frac{e^{-\lambda} \lambda^y}{y!}, \quad y = 0, 1, 2, \ldots\] Dalam bentuk eksponensial: \[f(y; \lambda) = \exp\left[y \log(\lambda) - \lambda - \log(y!)\right]\] di mana:

\(\theta = \log(\lambda)\). \(b(\theta) = e^\theta\). \(a(\phi) = 1\). \(c(y) = -\log(y!)\).

Parameter, Mean, dan Varians:

  • Parameter: \(\lambda\) (rata-rata kejadian, \(\lambda > 0\)).

  • Mean: \(E(Y) = \lambda\).

  • Varians: \(\text{Var}(Y) = \lambda\).

Distribusi Poisson digunakan dalam regresi Poisson untuk memodelkan frekuensi kejadian, seperti jumlah kasus penyakit dalam kelompok tertentu.

Model Regresi Logistik

Regresi logistik adalah model GLM yang digunakan untuk memodelkan variabel respons kategorik, terutama biner (sukses/gagal) atau multikategori, dengan menghubungkan probabilitas hasil ke prediktor melalui fungsi logit. Model ini cocok untuk data kategorik karena mengasumsikan respons mengikuti distribusi Bernoulli (untuk biner) atau Multinomial (untuk multikategori).

Untuk respons biner \(Y_i \in {0, 1}\) pada observasi ke-\(i\), regresi logistik memodelkan probabilitas sukses \(\pi_i = P(Y_i = 1)\) sebagai:

\[ \log\left(\frac{\pi_i}{1-\pi_i}\right) = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_p x_{ip} \] di mana:

\(\pi_i\): Probabilitas sukses untuk observasi ke-\(i\). \(\beta_0\): Intersep. \(\beta_j\): Koefisien untuk prediktor \(x_{ij}\) (\(j = 1, \ldots, p\)). \(x_{ij}\): Nilai prediktor ke-\(j\) untuk observasi ke-\(i\).

Fungsi link logit menghubungkan mean respons \(\mu_i = \pi_i\) ke prediktor linier \(\eta_i = \beta_0 + \sum \beta_j x_{ij}\), dan probabilitas dihitung sebagai:

\[ \pi_i = \frac{\exp(\eta_i)}{1 + \exp(\eta_i)} \] Untuk data multikategori dengan \(k\) kategori, regresi logistik multinomial memodelkan probabilitas kategori \(j\) relatif terhadap kategori referensi \(k\):

\[ \log\left(\frac{P(Y_i = j)}{P(Y_i = k)}\right) = \beta_{j0} + \beta_{j1} x_{i1} + \cdots + \beta_{jp} x_{ip}, \quad j = 1, \ldots, k-1 \] Probabilitas masing-masing kategori adalah:

\[ P(Y_i = j) = \frac{\exp(\beta_{j0} + \sum \beta_{jj} x_{ij})}{1 + \sum_{m=1}^{k-1} \exp(\beta_{m0} + \sum \beta_{mj} x_{ij})} \] Untuk data ordinal, model odds proporsional digunakan:

\[ \log\left(\frac{P(Y_i \leq j)}{P(Y_i > j)}\right) = \alpha_j + \beta_1 x_{i1} + \cdots + \beta_p x_{ip} \] di mana \(\alpha_j\) adalah ambang batas untuk kategori \(j\).

Regresi logistik digunakan untuk memprediksi probabilitas kejadian biner, seperti risiko penyakit berdasarkan faktor usia dan kebiasaan merokok, atau untuk mengklasifikasikan respons multikategori, seperti preferensi konsumen. Keunggulannya adalah kemampuan menangani prediktor campuran dan interpretasi koefisien sebagai log-odds.

Model Regresi Poisson

Regresi Poisson adalah model GLM yang digunakan untuk memodelkan variabel respons berupa hitungan atau frekuensi kejadian, yang diasumsikan mengikuti distribusi Poisson. Model ini cocok untuk data kategorik berupa jumlah kejadian dalam kategori, seperti jumlah kasus penyakit per wilayah.

Untuk respons \(Y_i\) (jumlah kejadian pada observasi ke-\(i\)), regresi Poisson memodelkan rata-rata kejadian \(\lambda_i = E(Y_i)\) dengan fungsi link log:

\[ \log(\lambda_i) = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_p x_{ip} \] di mana:

\(\lambda_i\): Rata-rata kejadian untuk observasi ke-\(i\). \(\beta_0\): Intersep. \(\beta_j\): Koefisien untuk prediktor \(x_{ij}\). \(x_{ij}\): Nilai prediktor ke-\(j\) untuk observasi ke-\(i\).

Fungsi link log memastikan \(\lambda_i > 0\), dan probabilitas distribusi Poisson adalah:

\[ P(Y_i = y_i) = \frac{e^{-\lambda_i} \lambda_i^{y_i}}{y_i!}, \quad y_i = 0, 1, 2, \ldots \]

Koefisien \(\beta_j\) diinterpretasikan sebagai perubahan log-rata-rata kejadian per unit perubahan prediktor \(x_j\), sehingga \(\exp(\beta_j)\) adalah rasio rata-rata kejadian.

Regresi Poisson digunakan untuk menganalisis frekuensi kejadian, seperti jumlah kecelakaan lalu lintas per wilayah berdasarkan faktor seperti kepadatan populasi, atau jumlah kasus penyakit berdasarkan faktor risiko. Keunggulannya adalah kemampuan menangani data hitungan dan interpretasi yang intuitif, tetapi model ini sensitif terhadap overdispersi (varians lebih besar dari mean).

Generalized Linear Model menyediakan kerangka fleksibel untuk analisis data kategorik melalui keluarga eksponensial, yang mencakup distribusi Bernoulli, Binomial, Multinomial, dan Poisson. Regresi logistik memodelkan probabilitas respons biner atau multikategori dengan fungsi link logit, cocok untuk klasifikasi dan prediksi. Regresi Poisson memodelkan frekuensi kejadian dengan fungsi link log, ideal untuk data hitungan. Pemilihan model tergantung pada jenis respons dan struktur data, memastikan analisis yang akurat dan relevan.

Inferensi pada Generalized Linear Model

Generalized Linear Model (GLM) adalah kerangka statistik untuk menganalisis data kategorik, seperti data biner atau hitungan, dengan menghubungkan mean respons ke prediktor melalui fungsi link. Inferensi pada GLM melibatkan estimasi parameter, evaluasi ekspektasi dan varians, serta diagnostik model untuk memastikan kesesuaian dengan data. Bagian ini menjelaskan lima aspek inferensi pada GLM: menduga ekspektasi dan varians, menaksir parameter, diagnostik model, serta metode estimasi dan inferensi untuk regresi logistik dan regresi Poisson. Semua perhitungan dilakukan secara manual dan dengan sintaks R tanpa package statistik, menggunakan data kecelakaan lalu lintas berikut:

Waktu Berkendara Jenis Kendaraan Cedera Tidak Cedera
Malam Motor 43 134
Malam Mobil 26 149
Siang Motor 29 23
Siang Mobil 22 26

Menduga Ekspektasi dan Varians dalam GLM

Ekspektasi (mean) dan varians dalam GLM ditentukan oleh distribusi respons dari keluarga eksponensial dan fungsi link.

Prosedur:

Ekspektasi: Mean respons \(\mu_i = E(Y_i)\) adalah probabilitas sukses \(\pi_i\) untuk data biner (Bernoulli) atau rata-rata kejadian \(\lambda_i\) untuk hitungan (Poisson).

Fungsi Link: Menghubungkan \(\mu_i\) ke prediktor linier \(\eta_i = \beta_0 + \sum \beta_j x_{ij}\), misalnya, logit untuk regresi logistik atau log untuk regresi Poisson.

Varians: Untuk Bernoulli, \(\text{Var}(Y_i) = \pi_i (1-\pi_i)\). Untuk Poisson, \(\text{Var}(Y_i) = \lambda_i\).

Perhitungan Manual: Untuk data kecelakaan, respons adalah cedera (1 = cedera, 0 = tidak cedera), mengikuti distribusi Bernoulli. Ekspektasi dan varians dihitung untuk setiap kelompok:

  • Malam, Motor:

Total: \(43 + 134 = 177\) \(\pi_1 = 43 / 177 \approx 0,24294\) Varians: \(\pi_1 (1 - \pi_1) = 0,24294 \times (1 - 0,24294) \approx 0,18397\)

  • Malam, Mobil:

Total: \(26 + 149 = 175\) \(\pi_2 = 26 / 175 \approx 0,14857\) Varians: \(0,14857 \times (1 - 0,14857) \approx 0,12649\)

  • Siang, Motor:

Total: \(29 + 23 = 52\) \(\pi_3 = 29 / 52 \approx 0,55769\) Varians: \(0,55769 \times (1 - 0,55769) \approx 0,24678\)

  • Siang, Mobil:

Total: \(22 + 26 = 48\) \(\pi_4 = 22 / 48 \approx 0,45833\) Varians: \(0,45833 \times (1 - 0,45833) \approx 0,24826\)

Untuk regresi Poisson, jumlah cedera per kelompok adalah \(\lambda_i\), misalnya, \(\lambda_1 = 43\) untuk malam-motor, dengan varians \(\text{Var}(Y_1) = 43\).

Sintaks R Manual:

# Data
data <- list(
  list(waktu="Malam", kendaraan="Motor", cedera=43, tidak_cedera=134),
  list(waktu="Malam", kendaraan="Mobil", cedera=26, tidak_cedera=149),
  list(waktu="Siang", kendaraan="Motor", cedera=29, tidak_cedera=23),
  list(waktu="Siang", kendaraan="Mobil", cedera=22, tidak_cedera=26)
)

# Menghitung ekspektasi dan varians
for (d in data) {
  total <- d$cedera + d$tidak_cedera
  pi <- d$cedera / total
  var <- pi * (1 - pi)
  cat(sprintf("%s, %s: pi = %.5f, var = %.5f\n", d$waktu, d$kendaraan, pi, var))
}
## Malam, Motor: pi = 0.24294, var = 0.18392
## Malam, Mobil: pi = 0.14857, var = 0.12650
## Siang, Motor: pi = 0.55769, var = 0.24667
## Siang, Mobil: pi = 0.45833, var = 0.24826
# Pengambilan keputusan
for (d in data) {
  pi <- d$cedera / (d$cedera + d$tidak_cedera)
  if (pi > 0.5) {
    cat(sprintf("%s, %s: Probabilitas cedera tinggi (%.3f), perhatikan keselamatan.\n", d$waktu, d$kendaraan, pi))
  } else {
    cat(sprintf("%s, %s: Probabilitas cedera rendah (%.3f), relatif aman.\n", d$waktu, d$kendaraan, pi))
  }
}
## Malam, Motor: Probabilitas cedera rendah (0.243), relatif aman.
## Malam, Mobil: Probabilitas cedera rendah (0.149), relatif aman.
## Siang, Motor: Probabilitas cedera tinggi (0.558), perhatikan keselamatan.
## Siang, Mobil: Probabilitas cedera rendah (0.458), relatif aman.

Menaksir Parameter

Parameter GLM, yaitu koefisien \(\beta\), diestimasi menggunakan maximum likelihood estimation (MLE) dengan memaksimalkan log-likelihood.

**Prosedur:*

  • Likelihood: Produk probabilitas observasi berdasarkan distribusi respons.

  • Log-Likelihood: Logaritma dari likelihood untuk mempermudah optimasi.

  • Optimasi: Menggunakan iterasi seperti Newton-Raphson untuk menemukan \(\hat{\beta}\).

  • Standar Error: Diperoleh dari invers matriks Hessian log-likelihood.

Perhitungan Manual (Regresi Logistik) Untuk regresi logistik, modelnya adalah:

\[ \log\left(\frac{\pi_i}{1-\pi_i}\right) = \beta_0 + \beta_1 \text{Waktu}{\text{Siang}} + \beta_2 \text{Kendaraan}{\text{Mobil}} \]

Data diubah ke format observasi individu:

Malam, Motor: 43 cedera (1), 134 tidak cedera (0). Malam, Mobil: 26 cedera (1), 149 tidak cedera (0). Siang, Motor: 29 cedera (1), 23 tidak cedera (0). Siang, Mobil: 22 cedera (1), 26 tidak cedera (0).

Total observasi: \(177 + 175 + 52 + 48 = 452\). Log-likelihood adalah:

\[ \ell(\beta) = \sum_{i=1}^n \left[ y_i \log(\pi_i) + (1-y_i) \log(1-\pi_i) \right] \]

dengan \(\pi_i = \exp(\eta_i) / (1 + \exp(\eta_i))\), dan \(\eta_i = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2}\), di mana \(x_{i1} = 1\) untuk siang, 0 untuk malam; \(x_{i2} = 1\) untuk mobil, 0 untuk motor.

Inisialisasi \(\beta_0 = -1, \beta_1 = 1, \beta_2 = -0,5\). Iterasi manual:

Hitung \(\eta_i\) dan \(\pi_i\).

Hitung gradien: \(\partial \ell / \partial \beta_j = \sum (y_i - \pi_i) x_{ij}\).

Perbarui \(\beta_j \leftarrow \beta_j + \alpha \cdot \text{gradien}\), dengan \(\alpha = 0,01\).

Setelah iterasi, misalkan konvergensi ke \(\hat{\beta}_0 \approx -0,876\), \(\hat{\beta}_1 \approx 1,234\), \(\hat{\beta}_2 \approx -0,623\).

Sintaks R Manual:

# Data individu
data_indiv <- c()
for (d in data) {
  for (i in 1:d$cedera) {
    data_indiv <- rbind(data_indiv, c(ifelse(d$waktu == "Siang", 1, 0), ifelse(d$kendaraan == "Mobil", 1, 0), 1))
  }
  for (i in 1:d$tidak_cedera) {
    data_indiv <- rbind(data_indiv, c(ifelse(d$waktu == "Siang", 1, 0), ifelse(d$kendaraan == "Mobil", 1, 0), 0))
  }
}
colnames(data_indiv) <- c("x1", "x2", "y")

# Fungsi sigmoid
sigmoid <- function(x) {
  1 / (1 + exp(-x))
}

# Log-likelihood
log_likelihood <- function(beta, data) {
  ll <- 0
  for (i in 1:nrow(data)) {
    eta <- beta[1] + beta[2] * data[i, "x1"] + beta[3] * data[i, "x2"]
    pi <- sigmoid(eta)
    ll <- ll + data[i, "y"] * log(pi + 1e-10) + (1 - data[i, "y"]) * log(1 - pi + 1e-10)
  }
  return(ll)
}

# Gradien
gradient <- function(beta, data) {
  grad <- c(0, 0, 0)
  for (i in 1:nrow(data)) {
    eta <- beta[1] + beta[2] * data[i, "x1"] + beta[3] * data[i, "x2"]
    pi <- sigmoid(eta)
    error <- data[i, "y"] - pi
    grad[1] <- grad[1] + error
    grad[2] <- grad[2] + error * data[i, "x1"]
    grad[3] <- grad[3] + error * data[i, "x2"]
  }
  return(grad)
}

# Iterasi MLE
beta <- c(-1, 1, -0.5)
learning_rate <- 0.01
for (i in 1:100) {
  grad <- gradient(beta, data_indiv)
  beta <- beta + learning_rate * grad
}

# Pengambilan keputusan
for (i in 2:3) {
  odds_ratio <- exp(beta[i])
  if (odds_ratio > 1) {
    cat(sprintf("Beta%d: Odds ratio = %.3f, meningkatkan risiko cedera.\n", i-1, odds_ratio))
  } else {
    cat(sprintf("Beta%d: Odds ratio = %.3f, menurunkan risiko cedera.\n", i-1, odds_ratio))
  }
}
## Beta1: Odds ratio = 4.327, meningkatkan risiko cedera.
## Beta2: Odds ratio = 0.581, menurunkan risiko cedera.

Diagnostik Model GLM

Diagnostik memastikan model GLM sesuai dengan data dan tidak melanggar asumsi.

Prosedur:

  • Residual Deviance: Mengukur kesesuaian model; nilai besar menunjukkan ketidaksesuaian.

  • Pearson Residual: \((y_i - \hat{\mu}_i) / \sqrt{\text{Var}(\hat{\mu}_i)}\) untuk mengidentifikasi observasi yang tidak sesuai.

  • Goodness-of-Fit: Uji seperti Hosmer-Lemeshow untuk regresi logistik.

  • Overdispersi: Untuk Poisson, memeriksa apakah varians sesuai dengan mean.

Untuk regresi logistik, Pearson residual dihitung per kelompok:

Malam, Motor: \(\hat{\pi}_1 = \exp(-0,876) / (1 + \exp(-0,876)) \approx 0,294\).

Residual: \((43/177 - 0,294) / \sqrt{0,294 \times 0,706 / 177} \approx -0,856\).

Malam, Mobil: \(\hat{\pi}_2 \approx 0,149\). Residual: \((26/175 - 0,149) / \sqrt{0,149 \times 0,851 / 175} \approx -0,027\).

Siang, Motor: \(\hat{\pi}_3 \approx 0,557\). Residual: \((29/52 - 0,557) / \sqrt{0,557 \times 0,443 / 52} \approx 0,032\).

Siang, Mobil: \(\hat{\pi}_4 \approx 0,294\). Residual: \((22/48 - 0,294) / \sqrt{0,294 \times 0,706 / 48} \approx 2,268\).

Sintaks R Manual:

# Pearson residual
for (d in data) {
  x1 <- ifelse(d$waktu == "Siang", 1, 0)
  x2 <- ifelse(d$kendaraan == "Mobil", 1, 0)
  eta <- beta[1] + beta[2] * x1 + beta[3] * x2
  pi <- sigmoid(eta)
  total <- d$cedera + d$tidak_cedera
  observed <- d$cedera / total
  residual <- (observed - pi) / sqrt(pi * (1 - pi) / total)
  if (abs(residual) > 2) {
    cat(sprintf("%s, %s: Residual besar (%.3f), model mungkin tidak fit.\n", d$waktu, d$kendaraan, residual))
  } else {
    cat(sprintf("%s, %s: Residual kecil (%.3f), model cukup fit.\n", d$waktu, d$kendaraan, residual))
  }
}
## Malam, Motor: Residual kecil (0.156), model cukup fit.
## Malam, Mobil: Residual kecil (-0.185), model cukup fit.
## Siang, Motor: Residual kecil (-0.248), model cukup fit.
## Siang, Mobil: Residual kecil (0.257), model cukup fit.

Metode Estimasi dan Inferensi Regresi Logistik

Estimasi: Regresi logistik memodelkan probabilitas cedera menggunakan distribusi Bernoulli dan fungsi link logit. MLE memaksimalkan log-likelihood dengan iterasi.

Inferensi:

Uji Wald: Statistik \(z = \hat{\beta}_j / \text{SE}(\hat{\beta}_j)\) untuk menguji signifikansi.

Odds Ratio: \(\exp(\hat{\beta}_j)\) menunjukkan perubahan odds per unit prediktor.

Likelihood Ratio Test: Membandingkan model penuh dan tereduksi.

Perhitungan Manual:

Odds Ratio: Untuk \(\hat{\beta}_1 \approx 1,234\), \(\exp(1,234) \approx 3,435\), berarti odds cedera di siang hari 3,435 kali lebih tinggi dibandingkan malam.

Standar Error: Diperoleh dari akar invers diagonal matriks Hessian, misalkan \(\text{SE}(\hat{\beta}_1) \approx 0,3\).

Uji Wald: \(z = 1,234 / 0,3 \approx 4,113\), signifikan (p-value < 0,05).

Sintaks R Manual:

# Odds ratio dan uji Wald
se <- c(0.3, 0.3, 0.3)  # Asumsi standar error
for (i in 1:3) {
  odds_ratio <- exp(beta[i])
  z <- beta[i] / se[i]
  p_value <- 2 * (1 - pnorm(abs(z)))
  if (p_value < 0.05) {
    cat(sprintf("Beta%d: Odds ratio = %.3f, signifikan (p = %.3f).\n", i-1, odds_ratio, p_value))
  } else {
    cat(sprintf("Beta%d: Odds ratio = %.3f, tidak signifikan (p = %.3f).\n", i-1, odds_ratio, p_value))
  }
}
## Beta0: Odds ratio = 0.312, signifikan (p = 0.000).
## Beta1: Odds ratio = 4.327, signifikan (p = 0.000).
## Beta2: Odds ratio = 0.581, tidak signifikan (p = 0.071).

Metode Estimasi dan Inferensi Regresi Poisson

Estimasi: Regresi Poisson memodelkan jumlah cedera menggunakan distribusi Poisson dan fungsi link log. MLE memaksimalkan log-likelihood.

Inferensi:

Uji Wald: Menguji signifikansi koefisien.

Rate Ratio: \(\exp(\hat{\beta}_j)\) menunjukkan perubahan rata-rata kejadian.

Likelihood Ratio Test: Mengevaluasi signifikansi prediktor.

Untuk regresi Poisson, modelnya adalah:

\[ \log(\lambda_i) = \beta_0 + \beta_1 \text{Waktu}{\text{Siang}} + \beta_2 \text{Kendaraan}{\text{Mobil}} \]

Log-likelihood Poisson:

\[ \ell(\beta) = \sum_{i=1}^n \left[ y_i \log(\lambda_i) - \lambda_i - \log(y_i!) \right] \]

dengan \(\lambda_i = \exp(\eta_i)\). Inisialisasi \(\beta_0 = 3,5\), \(\beta_1 = 0\), \(\beta_2 = 0\). Setelah iterasi, misalkan \(\hat{\beta}_0 \approx 3,76\), \(\hat{\beta}_1 \approx 0,21\), \(\hat{\beta}_2 \approx -0,43\).

Rate Ratio: \(\exp(0,21) \approx 1,234\), berarti rata-rata cedera di siang hari 1,234 kali lebih tinggi.

Uji Wald: Dengan \(\text{SE}(\hat{\beta}_1) \approx 0,15\), \(z = 0,21 / 0,15 \approx 1,4\).

Sintaks R Manual

# Regresi Poisson
poisson_log_likelihood <- function(beta, data) {
  ll <- 0
  for (d in data) {
    x1 <- ifelse(d$waktu == "Siang", 1, 0)
    x2 <- ifelse(d$kendaraan == "Mobil", 1, 0)
    lambda <- exp(beta[1] + beta[2] * x1 + beta[3] * x2)
    y <- d$cedera
    ll <- ll + y * log(lambda + 1e-10) - lambda
  }
  return(ll)
}

poisson_gradient <- function(beta, data) {
  grad <- c(0, 0, 0)
  for (d in data) {
    x1 <- ifelse(d$waktu == "Siang", 1, 0)
    x2 <- ifelse(d$kendaraan == "Mobil", 1, 0)
    lambda <- exp(beta[1] + beta[2] * x1 + beta[3] * x2)
    y <- d$cedera
    error <- y - lambda
    grad[1] <- grad[1] + error
    grad[2] <- grad[2] + error * x1
    grad[3] <- grad[3] + error * x2
  }
  return(grad)
}

# Iterasi MLE
beta_poisson <- c(3.5, 0, 0)
learning_rate <- 0.01
for (i in 1:100) {
  grad <- poisson_gradient(beta_poisson, data)
  beta_poisson <- beta_poisson + learning_rate * grad
}

# Pengambilan keputusan
for (i in 2:3) {
  rate_ratio <- exp(beta_poisson[i])
  if (rate_ratio > 1) {
    cat(sprintf("Beta%d: Rate ratio = %.3f, meningkatkan jumlah cedera.\n", i-1, rate_ratio))
  } else {
    cat(sprintf("Beta%d: Rate ratio = %.3f, menurunkan jumlah cedera.\n", i-1, rate_ratio))
  }
}
## Beta1: Rate ratio = 0.739, menurunkan jumlah cedera.
## Beta2: Rate ratio = 0.667, menurunkan jumlah cedera.

Inferensi pada GLM melibatkan estimasi ekspektasi dan varians berdasarkan distribusi respons, penaksiran parameter dengan MLE, diagnostik model untuk memastikan kesesuaian, serta metode estimasi dan inferensi spesifik untuk regresi logistik dan Poisson. Dalam kasus kecelakaan lalu lintas, regresi logistik memodelkan probabilitas cedera dengan distribusi Bernoulli, sedangkan regresi Poisson memodelkan jumlah cedera dengan distribusi Poisson. Perhitungan manual dan sintaks R tanpa package, dengan struktur if-else untuk interpretasi, memungkinkan pengambilan keputusan yang jelas berdasarkan odds ratio untuk logistik dan rate ratio untuk Poisson.

Regresi Logistik dengan Prediktor Nominal, Ordinal, dan Rasio

Regresi logistik adalah salah satu teknik statistik yang digunakan untuk menganalisis hubungan antara sebuah variabel respons biner (yang memiliki dua kategori) dengan satu atau lebih variabel prediktor. Dalam penggunaannya, variabel-variabel prediktor bisa memiliki berbagai jenis skala pengukuran, yaitu :

Dataset

Dataset yang digunakan adalah dataset Diabetes dengan 300 observasi. Variabel respon Diabetes Outcome merupakan variabel biner (0 = tidak, 1 = ya). Variabel prediktor terdiri dari:

  • Nominal : Sex
  • Ordinal: Social Class
  • Rasio: Blood Pressure, Glucose, Age

Dataset ini telah dibersihkan dan tidak mengandung nilai kosong.

Eksplorasi Data

Eksplorasi awal dilakukan dengan melihat struktur data dan ringkasan statistik:

Diabetes <- data.frame(
  Age = c(27, 47, 37, 43, 28, 39, 35, 34, 47, 48, 64, 35, 45, 60, 51, 38, 27, 26, 53, 54, 45, 52, 60, 58, 47, 54, 48, 47, 38, 59, 44, 42, 61, 40, 22, 33, 51, 25, 65, 50, 54, 74, 47, 61, 60, 51, 58, 44, 62, 60),
  Sex = factor(c("F", "F", "M", "M", "M", "F", "F", "M", "F", "F", "F", "M", "M", "F", "F", "F", "F", "F", "M", "M", "M", "F", "M", "F", "M", "M", "F", "F", "M", "M", "M", "M", "M", "M", "F", "M", "M", "F", "M", "M", "M", "F", "M", "F", "M", "F", "M", "M", "F", "M")),
  SocialClass = factor(c("MID", "LOW", "HIGH", "HIGH", "LOW", "LOW", "HIGH", "HIGH", "MID", "LOW", "MID", "MID", "MID", "LOW", "MID", "HIGH", "LOW", "LOW", "MID", "MID", "HIGH", "MID", "MID", "LOW", "MID", "LOW", "HIGH", "HIGH", "MID", "MID", "LOW", "LOW", "HIGH", "HIGH", "LOW", "MID", "MID", "HIGH", "LOW", "LOW", "MID", "LOW", "HIGH", "LOW", "LOW", "MID", "MID", "LOW", "MID", "HIGH")),
  BloodPressure = c(70, 96, 92, 80, 70, 80, 98, 63, 78, 98, 76, 88, 85, 78, 77, 64, 89, 94, 63, 59, 91, 58, 77, 76, 87, 87, 66, 90, 90, 81, 84, 88, 60, 64, 58, 73, 92, 90, 76, 60, 57, 97, 65, 68, 78, 80, 80, 94, 66, 87),
  Glucose = c(148, 85, 183, 114, 135, 110, 137, 78, 107, 133, 103, 174, 176, 92, 98, 134, 118, 112, 192, 174, 150, 164, 72, 84, 76, 173, 154, 143, 115, 165, 86, 135, 75, 128, 91, 98, 134, 110, 130, 174, 164, 98, 139, 135, 88, 70, 65, 136, 141, 166),
  DiabetesOutcome = factor(c(1, 0, 0, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 1, 1, 1, 0, 0, 0, 0, 1, 0, 0, 1, 1))
)
str(Diabetes)
## 'data.frame':    50 obs. of  6 variables:
##  $ Age            : num  27 47 37 43 28 39 35 34 47 48 ...
##  $ Sex            : Factor w/ 2 levels "F","M": 1 1 2 2 2 1 1 2 1 1 ...
##  $ SocialClass    : Factor w/ 3 levels "HIGH","LOW","MID": 3 2 1 1 2 2 1 1 3 2 ...
##  $ BloodPressure  : num  70 96 92 80 70 80 98 63 78 98 ...
##  $ Glucose        : num  148 85 183 114 135 110 137 78 107 133 ...
##  $ DiabetesOutcome: Factor w/ 2 levels "0","1": 2 1 1 1 2 1 2 1 1 2 ...
summary(Diabetes)
##       Age        Sex    SocialClass BloodPressure      Glucose     
##  Min.   :22.00   F:22   HIGH:13     Min.   :57.00   Min.   : 65.0  
##  1st Qu.:38.25   M:28   LOW :18     1st Qu.:66.50   1st Qu.: 98.0  
##  Median :47.00          MID :19     Median :79.00   Median :131.5  
##  Mean   :47.00                      Mean   :78.36   Mean   :125.2  
##  3rd Qu.:57.00                      3rd Qu.:88.75   3rd Qu.:149.5  
##  Max.   :74.00                      Max.   :98.00   Max.   :192.0  
##  DiabetesOutcome
##  0:25           
##  1:25           
##                 
##                 
##                 
##

Interpretasi :

  • Usia responden berkisar antara 22 hingga 74 tahun, dengan rata-rata dan median sebesar 47 tahun. Kuartil pertama (Q1) adalah 38.25 dan kuartil ketiga (Q3) adalah 57 tahun.
  • Jenis kelamin terdiri dari 22 perempuan (F) dan 28 laki-laki (M).
  • Kelas sosial terbagi menjadi 13 orang dari kelas tinggi (HIGH), 19 dari kelas menengah (MID), dan 18 dari kelas rendah (LOW).
  • Tekanan darah berkisar antara 57 hingga 98, dengan rata-rata 78.36 dan median 79.
  • Kadar glukosa darah berada dalam rentang 65 hingga 192, dengan rata-rata 125.2 dan median 131.5.
  • Hasil diagnosis diabetes (DiabetesOutcome) menunjukkan jumlah yang seimbang, yaitu 25 orang tidak menderita diabetes (0) dan 25 orang menderita diabetes (1).

Model Regresi Logistik

Model regresi logistik dibentuk sebagai berikut:

\[ \log\left(\frac{p}{1 - p}\right) = \beta_0 + \beta_1 \cdot \text{Age} + \beta_2 \cdot \text{Sex} + \beta_3 \cdot \text{SocialClass} + \beta_4 \cdot \text{BloddPressure} + \beta_5 \cdot \text{Glucose} \]

Model dapat diestimasi dengan perintah berikut :

model_log <- glm(DiabetesOutcome ~ Age + + Sex + SocialClass + BloodPressure + Glucose,
                 data = Diabetes, family = binomial)
summary(model_log)
## 
## Call:
## glm(formula = DiabetesOutcome ~ Age + +Sex + SocialClass + BloodPressure + 
##     Glucose, family = binomial, data = Diabetes)
## 
## Coefficients:
##                 Estimate Std. Error z value Pr(>|z|)  
## (Intercept)     0.576783   2.843592   0.203   0.8393  
## Age             0.003693   0.026688   0.138   0.8900  
## SexM           -1.687655   0.730564  -2.310   0.0209 *
## SocialClassLOW -0.311091   0.831728  -0.374   0.7084  
## SocialClassMID -0.079371   0.841407  -0.094   0.9248  
## BloodPressure  -0.033862   0.027108  -1.249   0.2116  
## Glucose         0.023900   0.010611   2.252   0.0243 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 69.315  on 49  degrees of freedom
## Residual deviance: 57.812  on 43  degrees of freedom
## AIC: 71.812
## 
## Number of Fisher Scoring iterations: 4

Interpretasi :

  • Intercept = 0.577
    Nilai intercept tidak signifikan (p = 0.839), artinya probabilitas dasar (saat semua prediktor = 0) tidak berbeda signifikan dari 0.5 (dalam konteks logit scale).
  • Age = 0.0037 (p = 0.890)
    Usia tidak signifikan dalam memprediksi diabetes. Setiap kenaikan 1 tahun usia tidak memberikan pengaruh berarti terhadap kemungkinan menderita diabetes.
  • SexM = -1.688 (p = 0.0209)
    Jenis kelamin laki-laki (M) secara signifikan mengurangi peluang terkena diabetes dibandingkan perempuan (referensi: Female).
    Odds laki-laki untuk terkena diabetes adalah sekitar exp(-1.688) ≈ 0.19 kali dibandingkan perempuan.

  • SocialClassLOW & MID (vs referensi HIGH)

    • LOW = -0.311 (p = 0.708)

    • MID = -0.079 (p = 0.925)
      Kelas sosial tidak signifikan dalam model ini. Tidak ada perbedaan yang berarti antara kelas sosial terhadap risiko diabetes.

  • BloodPressure = -0.034 (p = 0.211)
    Tekanan darah tidak signifikan memengaruhi risiko diabetes.
  • Glucose = 0.0239 (p = 0.0243)
    Kadar glukosa secara signifikan meningkatkan peluang terkena diabetes. Setiap kenaikan 1 unit glukosa meningkatkan odds terkena diabetes sebesar exp(0.0239) ≈ 1.024 atau naik sekitar 2.4%.

Evaluasi Model

Kurva 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, var
roc_curve <- roc(Diabetes$DiabetesOutcome, fitted(model_log))
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
plot(roc_curve, main = "Kurva ROC")

auc(roc_curve)
## Area under the curve: 0.7424

Nilai 0.7424 berarti bahwa model memiliki kemampuan yang cukup baik dalam membedakan antara individu yang menderita diabetes dan yang tidak.

Uji Hosmer-Lemeshow

library(ResourceSelection)
## Warning: package 'ResourceSelection' was built under R version 4.3.3
## ResourceSelection 0.3-6   2023-06-27
hoslem.test(Diabetes$DiabetesOutcome, fitted(model_log))
## Warning in Ops.factor(1, y): '-' not meaningful for factors
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  Diabetes$DiabetesOutcome, fitted(model_log)
## X-squared = 50, df = 8, p-value = 4.087e-08

Hasil evaluasi menunjukkan bahwa model tidak memiliki ketepatan klasifikasi yang baik dan terdapat masalah kecocokan model berdasarkan uji Hosmer-Lemeshow (p-value < 0.05).

Pemilihan Model Regresi Logistik dan Evaluasi

Membangun Model Regresi Logistik: Pendekatan Confirmatory dan Exploratory

Dalam analisis regresi logistik, pemilihan pendekatan dalam membangun model sangat penting untuk memperoleh model yang akurat dalam memperkirakan probabilitas kejadian suatu outcome biner. Terdapat dua pendekatan utama, yaitu Confirmatory dan Exploratory.

Confirmatory

Pendekatan ini digunakan ketika peneliti sudah memiliki dugaan atau teori sebelumnya mengenai hubungan antara variabel-variabel prediktor dan variabel respons. Fokus utamanya adalah menguji hipotesis yang sudah dirumuskan, bukan sekadar mencari model dengan performa terbaik.

Karakteristik :

  • Model dikembangkan berdasarkan teori ilmiah atau hasil studi terdahulu.

  • Tujuannya adalah mengonfirmasi apakah pengaruh suatu variabel benar-benar signifikan.

  • Peneliti biasanya membangun model penuh terlebih dahulu, lalu menilai apakah menambah atau menghapus variabel (termasuk interaksi) meningkatkan performa model secara signifikan.

  • Uji statistik seperti Likelihood Ratio Test (LRT) digunakan untuk membandingkan model yang berbeda.

Contoh :
Seorang peneliti ingin mengetahui apakah tingkat stres kerja dan durasi tidur memengaruhi kemungkinan seseorang mengalami tekanan darah tinggi. Berdasarkan teori kesehatan, model langsung dibangun menggunakan kedua variabel tersebut. Selanjutnya diuji apakah variabel “durasi tidur” secara signifikan berkontribusi dalam model.

Exploratory

Pendekatan ini dipakai ketika peneliti belum memiliki teori yang pasti, dan ingin menggali kemungkinan hubungan yang terdapat dalam data. Tujuan utamanya adalah menemukan model dengan kombinasi prediktor terbaik secara statistik.

Karakteristik :

  • Proses pemodelan bersifat bertahap dan fleksibel.

  • Pemilihan variabel dilakukan berdasarkan metrik statistik seperti AIC, deviance, atau log-likelihood.

  • Tiga teknik umum digunakan:

    • Forward Selection: Dimulai dari model kosong, variabel ditambahkan satu per satu.

    • Backward Elimination: Dimulai dari model penuh, variabel yang tidak signifikan dihapus.

    • Stepwise Selection: Kombinasi forward dan backward, memungkinkan variabel masuk dan keluar model secara dinamis.

Tujuan :
Mendapatkan model yang ringkas namun efektif (parsimonious) dalam memprediksi outcome.

Contoh :
Seorang analis data memiliki dataset pelanggan dengan variabel seperti usia, pengeluaran, frekuensi belanja, dan jenis kelamin, dan ingin memprediksi apakah pelanggan akan berhenti berlangganan. Karena tidak ada teori baku, ia menggunakan stepwise selection untuk menemukan kombinasi variabel terbaik yang mampu memprediksi churn pelanggan secara akurat.

Seleksi Model dengan Stepwise

Langkah ini menggunakan model awal regresi logistik, kemudian dilakukan seleksi otomatis dengan pendekatan dua arah (forward dan backward) :

step_model <- step(model_log, direction = "both")
## Start:  AIC=71.81
## DiabetesOutcome ~ Age + +Sex + SocialClass + BloodPressure + 
##     Glucose
## 
##                 Df Deviance    AIC
## - SocialClass    2   57.973 67.973
## - Age            1   57.831 69.831
## - BloodPressure  1   59.431 71.431
## <none>               57.812 71.812
## - Glucose        1   63.776 75.776
## - Sex            1   63.991 75.991
## 
## Step:  AIC=67.97
## DiabetesOutcome ~ Age + Sex + BloodPressure + Glucose
## 
##                 Df Deviance    AIC
## - Age            1   57.985 65.985
## - BloodPressure  1   59.778 67.778
## <none>               57.973 67.973
## + SocialClass    2   57.812 71.812
## - Sex            1   63.993 71.993
## - Glucose        1   64.232 72.232
## 
## Step:  AIC=65.98
## DiabetesOutcome ~ Sex + BloodPressure + Glucose
## 
##                 Df Deviance    AIC
## - BloodPressure  1   59.818 65.818
## <none>               57.985 65.985
## + Age            1   57.973 67.973
## + SocialClass    2   57.831 69.831
## - Sex            1   64.006 70.006
## - Glucose        1   64.240 70.240
## 
## Step:  AIC=65.82
## DiabetesOutcome ~ Sex + Glucose
## 
##                 Df Deviance    AIC
## <none>               59.818 65.818
## + BloodPressure  1   57.985 65.985
## + Age            1   59.778 67.778
## - Sex            1   65.321 69.321
## + SocialClass    2   59.472 69.472
## - Glucose        1   66.362 70.362
summary(step_model)
## 
## Call:
## glm(formula = DiabetesOutcome ~ Sex + Glucose, family = binomial, 
##     data = Diabetes)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)  
## (Intercept) -2.24042    1.25951  -1.779   0.0753 .
## SexM        -1.52534    0.69098  -2.208   0.0273 *
## Glucose      0.02443    0.01041   2.347   0.0189 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 69.315  on 49  degrees of freedom
## Residual deviance: 59.818  on 47  degrees of freedom
## AIC: 65.818
## 
## Number of Fisher Scoring iterations: 4

Berdasarkan hasil stepwise tersebut didapatkan bahwa jenis kelamin dan kadar glukosa merupakan prediktor signifikan terhadap outcome diabetes. Pria cenderung memiliki risiko diabetes yang lebih rendah dibanding wanita, pada kadar glukosa yang sama. Dapat disimpulkan bahwa model ini sudah lebih baik dari model kosong, namun deviance masih cukup besar, sehingga model bisa ditingkatkan lebih lanjut.

Evaluasi Model Akhir

Pseudo R-Squared

null_model <- glm(DiabetesOutcome ~ 1, data = Diabetes, family = binomial)
1 - (logLik(step_model) / logLik(null_model))
## 'log Lik.' 0.1370134 (df=3)

Nilai log-likelihood model sebesar 0.137 menunjukkan tingkat kecocokan model terhadap data yang telah diestimasi. Semakin tinggi (lebih mendekati nol atau positif) nilai log-likelihood, semakin baik model menjelaskan data dibandingkan model nol (tanpa prediktor).

Kurva ROC Model Akhir

roc_step <- roc(Diabetes$DiabetesOutcome, fitted(step_model))
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
plot(roc_step, main = "ROC Curve - Model Akhir")

auc(roc_step)
## Area under the curve: 0.7176

Kurva ROC menunjukkan bahwa model akhir memiliki performa yang cukup baik dalam membedakan antara penderita dan non-penderita diabetes. Kurva yang melengkung ke kiri atas menandakan kemampuan klasifikasi yang lebih baik daripada model acak, sehingga model layak digunakan dalam prediksi.

Distribusi Multinomial

Distribusi multinomial adalah perluasan dari distribusi binomial untuk lebih dari dua kategori. Jika \(X_1, X_2, \ldots, X_k\) menyatakan banyaknya kejadian dalam masing-masing dari \(k\) kategori, maka:

\[ P(X_1 = x_1, \ldots, X_k = x_k) = \frac{n!}{x_1! x_2! \ldots x_k!} p_1^{x_1} p_2^{x_2} \ldots p_k^{x_k} \]

dengan \(\sum_{i=1}^k x_i = n\) dan \(\sum_{i=1}^k p_i = 1\).

Studi Kasus

Sebuah survei dilakukan terhadap 10 orang yang diminta memilih salah satu dari tiga jenis motor favorit: Honda (A), Ducati (D), dan Yamaha (Y).

Hasil survei:

  • Honda : 4 orang
  • Ducati : 3 orang
  • Yamaha : 3 orang

Probabilitas teoritik preferensi:

  • \(p_H = 0.4\),
  • \(p_D = 0.3\),
  • \(p_Y = 0.3\)

Pertanyaannya: Berapa peluang bahwa dalam 20 orang akan ada 6 yang memilih Honda, 8 memilih Ducati, dan 6 memilih Yamaha?

Rumus Distribusi Multinomial Distribusi peluang multinomial:

\[ P(X_1 = x_1, \ldots, X_k = x_k) = \frac{n!}{x_1! x_2! \ldots x_k!} p_1^{x_1} p_2^{x_2} \ldots p_k^{x_k} \]

dengan:

  • \(n = 20\), \(x_1 = 6\), \(x_2 = 8\), \(x_3 = 6\)
  • \(p_1 = 0.4\), \(p_2 = 0.3\), \(p_3 = 0.3\)

Perhitungan Manual di R

n <- 20 
x <- c(6, 8, 6) 
p <- c(0.3, 0.4, 0.3)
faktorial_total <- factorial(n) 
faktorial_x <- prod(factorial(x))
koefisien <- faktorial_total / faktorial_x
peluang <- koefisien * prod(p^x);peluang
## [1] 0.0405391

Probabilitas bahwa 20 orang memilih 6 Honda, 8 Ducati, dan 6 Yamaha (dengan proporsi preferensi 0.3, 0.4, dan 0.3) adalah 0.0405391. Distribusi multinomial digunakan untuk menghitung peluang dalam percobaan dengan beberapa kategori hasil. Rumus dasarnya merupakan generalisasi dari binomial untuk lebih dari dua kategori.

Multinomial Logistik Regression

Model ini digunakan untuk memodelkan hubungan antara satu variabel respon kategorik (>2 kategori) dan satu atau lebih variabel prediktor.

Misalkan \(Y\) memiliki \(K\) kategori, dan kita pilih referensi (baseline) kategori \(K\), maka model logit untuk kategori \(j\) adalah:

\[ \log \left( \frac{P(Y = j)}{P(Y = K)} \right) = \beta_{j0} + \beta_{j1}x_1 + \ldots + \beta_{jp}x_p \]

untuk \(j = 1, 2, \ldots, K - 1\).

Baseline Category Logit Model

Baseline-category logit model adalah model regresi logistik untuk variabel respon kategorik dengan lebih dari dua kategori (nominal). Model ini menggunakan satu kategori sebagai acuan (baseline) dan membandingkan kategori lainnya terhadap baseline tersebut dalam bentuk logit.

\[ \log \left(\frac{\pi_j}{\pi_c}\right), \quad j=1,\ldots,c-1 \]

dengan:

  • \(\pi_j\) adalah probabilitas respon berada di kategori \(j\)
  • \(\pi_c\) adalah probabilitas respon berada di kategori acuan (baseline)

Maka, terdapat sebanyak \((c-1)\) fungsi logit.

Catatan: Kategori baseline bisa ditentukan secara eksplisit, tetapi default di R adalah kategori terakhir.

Model Regresi Jika terdapat satu prediktor \(x\), maka bentuk umum model logit-nya adalah:

\[ \log \left(\frac{\pi_j}{\pi_c}\right) = \alpha_j + \beta_j x, \quad j = 1, \ldots, c-1 \]

Contoh Kasus: 3 Kategori Respon Misalkan respon \(Y\) memiliki tiga kategori: \(Y \in \{1,2,3\}\), dan kita gunakan kategori ke-3 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 \]

Terdapat dua model logit, satu untuk perbandingan kategori 1 dengan 3, dan satu lagi untuk kategori 2 dengan 3.

Relasi Antar Kategori Jika ingin menghitung logit antara kategori 1 dan 2:

\[ \log \left(\frac{\pi_1}{\pi_2}\right) = \log \left(\pi_1 / \pi_3\right) - \log \left(\pi_2 / \pi_3\right) \]

\[ = (\alpha_1 + \beta_1 x) - (\alpha_2 + \beta_2 x) = (\alpha_1 - \alpha_2) + (\beta_1 - \beta_2)x \]

Model baseline-category logit:

  • Digunakan untuk respon dengan kategori \(> 2\)
  • Menghasilkan \((c-1)\) fungsi logit terhadap satu baseline
  • Logit antara kategori selain baseline dapat dihitung dari selisih dua logit terhadap baseline

Implementasi di R menggunakan fungsi multinom() dari package nnet, dan kategori baseline bisa diatur dengan relevel().

Estimasi Parameter

Estimasi dilakukan menggunakan metode maximum likelihood dengan algoritma iteratif seperti Newton-Raphson.

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)\) dan \(y_{ij} = 1\) jika \(Y_i = j\).

Contoh Kasus

Sebuah perusahaan pengembang aplikasi pembelajaran bahasa asing ingin mengetahui faktor-faktor yang memengaruhi jenis langganan yang dipilih pengguna setelah masa uji coba gratis selama 14 hari. Perusahaan ingin mengoptimalkan strategi pemasaran mereka berdasarkan kebiasaan pengguna selama masa uji coba.

Mereka melakukan survei terhadap 300 pengguna yang telah menyelesaikan masa uji coba, dan mengumpulkan data sebagai berikut :

  • Keputusan Langganan : Jenis paket berlangganan yang dipilih setelah masa uji coba (tidak berlangganan, paket standar, paket pro)

  • Durasi Penggunaan Aplikasi : Total jam belajar selama masa uji coba ( jam belajar aktif).

  • Bahasa Favorit :Bahasa yang paling sering dipelajari selama masa uji coba (Inggris, Jepang, Prancis).

  • Indeks Interaksi Fitur : Skor aktivitas pengguna dengan fitur aplikasi tambahan (latihan interaktif, forum diskusi, gamifikasi kuis).

Tujannya adalah mengetahui bagaimana durasi penggunaan aplikasi, bahasa favorit, dan interaksi fitur memengaruhi pilihan keputusan langganan pelanggan.

Simulasi Data

set.seed(55)
n <- 300  

# Variabel-variabel simulasi
FavoriteLanguage <- sample(c("Inggris", "Jepang", "Prancis"), n, replace = TRUE)

Duration <- round(rnorm(n, mean = 15, sd = 5))  # Total jam belajar selama masa uji coba
Duration <- pmax(Duration, 0)  # Tidak ada durasi negatif

FeatureInteraction <- round(pmax(rnorm(n, mean = 70, sd = 15), 0))  # Indeks interaksi fitur

# Simulasikan jenis langganan berdasarkan bahasa favorit
Subscription <- sapply(FavoriteLanguage, function(lang) {
  if (lang == "Inggris") {
    sample(c("Tidak Berlangganan", "Paket Standar", "Paket Pro"), size = 1,
           prob = c(0.3, 0.4, 0.3))
  } else if (lang == "Jepang") {
    sample(c("Tidak Berlangganan", "Paket Standar", "Paket Pro"), size = 1,
           prob = c(0.4, 0.3, 0.3))
  } else {
    sample(c("Tidak Berlangganan", "Paket Standar", "Paket Pro"), size = 1,
           prob = c(0.2, 0.3, 0.5))
  }
})

# Gabungkan menjadi data frame
df <- data.frame(
  Subscription = factor(Subscription),
  Duration = Duration,
  FavoriteLanguage = factor(FavoriteLanguage),
  FeatureInteraction = FeatureInteraction
)

# Set baseline untuk analisis regresi (misalnya: Tidak Berlangganan)
df$Subscription <- relevel(df$Subscription, ref = "Tidak Berlangganan")

# Tampilkan data awal
head(df)

Estimasi Model

# Load library
library(nnet)

# Estimasi model
model_multinom <- multinom(Subscription ~ Duration + FavoriteLanguage + FeatureInteraction, data = df)
## # weights:  18 (10 variable)
## initial  value 329.583687 
## iter  10 value 312.706609
## final  value 312.420301 
## converged
# Tampilkan ringkasan hasil estimasi
summary(model_multinom)
## Call:
## multinom(formula = Subscription ~ Duration + FavoriteLanguage + 
##     FeatureInteraction, data = df)
## 
## Coefficients:
##               (Intercept)     Duration FavoriteLanguageJepang
## Paket Pro       2.2988036 -0.037814635             -0.5249742
## Paket Standar   0.6973863  0.002894146             -0.5824085
##               FavoriteLanguagePrancis FeatureInteraction
## Paket Pro                   0.9807760       -0.022726521
## Paket Standar               0.2208183       -0.004250402
## 
## Std. Errors:
##               (Intercept)   Duration FavoriteLanguageJepang
## Paket Pro       0.9065429 0.03116649              0.3751310
## Paket Standar   0.8956554 0.03085098              0.3505482
##               FavoriteLanguagePrancis FeatureInteraction
## Paket Pro                   0.3997707         0.01038147
## Paket Standar               0.4008121         0.01001261
## 
## Residual Deviance: 624.8406 
## AIC: 644.8406

Nilai P-Value dan Interpretasi

# Menghitung p-value
z_val <- summary(model_multinom)$coefficients / summary(model_multinom)$standard.errors
p_val <- 2 * (1 - pnorm(abs(z_val)))

# Tampilkan koefisien, z-value, dan p-value
cbind(
  Coefficients = round(summary(model_multinom)$coefficients, 4),
  `Z value` = round(z_val, 3),
  `P value` = round(p_val, 4)
)
##               (Intercept) Duration FavoriteLanguageJepang
## Paket Pro          2.2988  -0.0378                -0.5250
## Paket Standar      0.6974   0.0029                -0.5824
##               FavoriteLanguagePrancis FeatureInteraction (Intercept) Duration
## Paket Pro                      0.9808            -0.0227       2.536   -1.213
## Paket Standar                  0.2208            -0.0043       0.779    0.094
##               FavoriteLanguageJepang FavoriteLanguagePrancis FeatureInteraction
## Paket Pro                     -1.399                   2.453             -2.189
## Paket Standar                 -1.661                   0.551             -0.425
##               (Intercept) Duration FavoriteLanguageJepang
## Paket Pro          0.0112   0.2250                 0.1617
## Paket Standar      0.4362   0.9253                 0.0966
##               FavoriteLanguagePrancis FeatureInteraction
## Paket Pro                      0.0142             0.0286
## Paket Standar                  0.5817             0.6712

Kesimpulan

Berdasarkan hasil regresi logistik multinomial, ditemukan beberapa faktor yang secara signifikan memengaruhi kecenderungan seseorang memilih paket langganan dibandingkan dengan kategori referensi.

  • Bahasa favorit Jepang berpengaruh negatif signifikan terhadap pemilihan Paket Pro (z = -1.399, p = 0.1617) dan Paket Standar (z = -1.661, p = 0.0966). Meskipun nilai p untuk Paket Standar mendekati batas signifikansi 10%, secara umum ini menunjukkan bahwa pengguna dengan bahasa Jepang cenderung kurang memilih kedua paket tersebut, terutama Paket Standar.

  • Sebaliknya, bahasa favorit Prancis menunjukkan pengaruh positif dan signifikan terhadap pemilihan Paket Pro (z = 2.453, p = 0.0142), yang berarti pengguna yang menyukai bahasa Prancis lebih mungkin memilih Paket Pro dibandingkan kategori referensi. Namun, pengaruh terhadap Paket Standar tidak signifikan (p = 0.5817).

  • Durasi penggunaan menunjukkan hubungan negatif signifikan terhadap pemilihan Paket Pro (z = -1.213, p = 0.2250), tetapi tidak signifikan secara statistik. Pengaruh terhadap Paket Standar juga tidak signifikan (p = 0.9253).

  • Interaksi fitur (FeatureInteraction) tidak menunjukkan pengaruh signifikan terhadap pemilihan baik Paket Pro (p = 0.0286) maupun Paket Standar (p = 0.6712). Meski p untuk Paket Pro di bawah 0.05, nilai z = -2.189 menunjukkan pengaruh negatif; namun sebaiknya diuji lebih lanjut untuk memastikan stabilitas efek ini.

Regresi Logistik Ordinal

Regresi logistik ordinal digunakan ketika variabel respon Y bersifat ordinal (memiliki urutan), misalnya tingkat kepuasan: Rendah, Sedang, Tinggi. Model ini berbeda dengan : • Regresi logistik biner: hanya 2 kategori • Regresi logistik multinomial: kategori > 2 tetapi tidak berurutan

Konsep Cumulative Logit Model

Model yang digunakan adalah Cumulative Logit Model dengan asumsi proportional odds :

\[ \log \left( \frac{P(Y \leq j)}{P(Y > j)} \right) = \alpha_j + \beta x \]

  • \(\alpha_j\): intercept khusus untuk kategori ke-\(j\)
  • \(\beta\): koefisien regresi (sama untuk semua kategori kumulatif)

Untuk \(c\) kategori, terdapat \((c - 1)\) model logit kumulatif.

Interpretasi Koefisien

Koefisien \(\beta\) menjelaskan efek \(x\) terhadap kemungkinan berada pada kategori yang lebih rendah atau sama.

  • Jika \(\beta > 0\): semakin besar \(x\), semakin tinggi peluang berada di kategori rendah.
  • Jika \(\beta < 0\): semakin besar \(x\), semakin besar peluang berada di kategori tinggi.

Odds ratio:

\[ \text{OR} = e^{\beta} \]

Studi Kasus

Misal kita memiliki data terkait tingkat kepuasan pelanggan (1 : Tidak Puas, 2: Cukup, 3: Puas) terhadap kualitsa pembelajaran :

set.seed(55)
n <- 200
speed <- round(runif(n, 1, 10))
satisfaction <- cut(5 + 0.5*speed + rnorm(n),
breaks = c(-Inf, 5.5, 7.5, Inf),
labels = c("Tidak Puas", "Cukup", "Puas"),
ordered_result = TRUE)
df <- data.frame(satisfaction, speed)
head(df)

Estimasi Model Ordinal

library(MASS)
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
model_ord <- polr(satisfaction ~ speed, data = df, Hess = TRUE)
summary(model_ord)
## Call:
## polr(formula = satisfaction ~ speed, data = df, Hess = TRUE)
## 
## Coefficients:
##        Value Std. Error t value
## speed 0.7743    0.09293   8.332
## 
## Intercepts:
##                  Value  Std. Error t value
## Tidak Puas|Cukup 0.7475 0.3883     1.9250 
## Cukup|Puas       3.7939 0.4964     7.6432 
## 
## Residual Deviance: 245.0161 
## AIC: 251.0161

Berdasarkan hasil regresi logistik ordinal terhadap variabel satisfaction dengan prediktor speed, diperoleh bahwa kecepatan (speed) berpengaruh positif dan signifikan terhadap tingkat kepuasan pengguna.

Koefisien untuk speed adalah 0.7743 dengan nilai t = 8.332, yang menunjukkan bahwa setiap kenaikan satu unit pada kecepatan akan meningkatkan peluang responden berada pada kategori kepuasan yang lebih tinggi (misalnya, dari Tidak Puas ke Cukup, atau dari Cukup ke Puas). Karena nilai t-nya besar dan tidak mendekati nol, ini menunjukkan bahwa pengaruh speed signifikan secara statistik.

Model ini juga menghasilkan dua ambang batas (cutpoints/intercepts):

“Tidak Puas|Cukup” sebesar 0.7475

“Cukup|Puas” sebesar 3.7939

Kedua nilai ini digunakan untuk menghitung probabilitas kumulatif kategori kepuasan. Misalnya, semakin besar nilai speed, semakin kecil kemungkinan responden berada di kategori Tidak Puas atau Cukup, dan semakin besar kemungkinan berada di kategori Puas.

Nilai Residual Deviance sebesar 245.0161 dan AIC sebesar 251.0161 memberikan informasi tentang kelayakan model, yang bisa dibandingkan dengan model lain jika diperlukan.

Nilai P-Value

(ctable <- coef(summary(model_ord)))
##                      Value Std. Error  t value
## speed            0.7742641  0.0929311 8.331593
## Tidak Puas|Cukup 0.7475341  0.3883392 1.924951
## Cukup|Puas       3.7938935  0.4963761 7.643183
p <- pnorm(abs(ctable[, "t value"]), lower.tail = FALSE) * 2
(ctable <- cbind(ctable, "p value" = round(p, 4)))
##                      Value Std. Error  t value p value
## speed            0.7742641  0.0929311 8.331593  0.0000
## Tidak Puas|Cukup 0.7475341  0.3883392 1.924951  0.0542
## Cukup|Puas       3.7938935  0.4963761 7.643183  0.0000

Karena p-value jauh lebih kecil dari 0.05, maka dapat disimpulkan bahwa speed berpengaruh secara signifikan terhadap kepuasan (satisfaction). Artinya, terdapat bukti statistik yang kuat bahwa peningkatan kecepatan layanan secara nyata meningkatkan kemungkinan responden merasa lebih puas.

Prediksi Probabilitas

newdata <- data.frame(speed = 5:9)
predict(model_ord, newdata = newdata, type = "probs")
##    Tidak Puas      Cukup      Puas
## 1 0.042136688 0.43851623 0.5193471
## 2 0.019878255 0.27919998 0.7009218
## 3 0.009263979 0.15512124 0.8356148
## 4 0.004292524 0.07886346 0.9168440
## 5 0.001983627 0.03815371 0.9598627

Goodness-of-Fit dan Proportional Odds

Model cumulative logit mengasumsikan efek prediktor sama untuk setiap cutof. Jika tidak, pertimbangkan model non-proportional odds seperti generalized ordinal model.

Alternatif Model Ordinal

Selain cumulative logit, model ordinal lainnya :

  • Adjacent-category logit

  • Continuation-ratio (sequential) logit

    Model tersebut dapat digunakan saat asumsi proportional odds tidak terpenuhi

Kesimpulan

  • Regresi ordinal efektif untuk respon berurutan.

  • Model cumulative logit menginterpretasikan efek dalam bentuk log-odds kumulatif.

  • Implementasi di R dengan fungsi polr() dari package MASS.

    Jika diperlukan validasi lanjut, bisa digunakan uji devians atau likelihood ratio test

Asumsi Paralelisme dalam Regresi Logistik Ordinal

Model regresi logistik ordinal yang paling umum digunakan adalah Cumulative Logit Model dengan asumsi Proportional Odds.

Asumsi ini dikenal juga sebagai asumsi paralelisme (parallel lines assumption).

Definisi Asumsi Paralelisme:
Asumsi paralelisme menyatakan bahwa koefisien regresi (β) sama untuk setiap kategori kumulatif dari variabel respon.

Bentuk umum model:

\[ \log \left( \frac{P(Y \leq j)}{P(Y > j)} \right) = \alpha_j + \beta x \]

Untuk \(j = 1, \ldots, c - 1\):

  • Hanya intercept (\(\alpha_j\)) yang berbeda-beda.
  • Koefisien \(\beta\) tetap sama untuk semua fungsi logit kumulatif.

Visualisasi Dalam asumsi paralelisme, kurva logit kumulatif dari tiap kategori terhadap prediktor akan memiliki kemiringan yang sama (paralel), hanya berbeda posisi (intercept).

Konsekuensi Pelanggaran Asumsi

Jika asumsi ini tidak terpenuhi :

  • Efek prediktor berbeda untuk setiap batas kategori.

  • Model cumulative logit tidak valid.

  • Perlu gunakan model alternatif: – Generalized Ordinal Logistic Regression – Partial Proportional Odds Model

    Pengujian Asumsi Paralelisme Untuk memeriksa validitas asumsi, dapat digunakan :

  • Likelihood Ratio Test antara model proportional dan non-proportional

  • Brant Test (paket brant di R)

library(MASS)
library(brant)
## Warning: package 'brant' was built under R version 4.3.3
model <- polr(satisfaction ~ speed, data = data, Hess = TRUE)
brant(model)
## -------------------------------------------- 
## Test for X2  df  probability 
## -------------------------------------------- 
## Omnibus      2.99    1   0.08
## speed        2.99    1   0.08
## -------------------------------------------- 
## 
## H0: Parallel Regression Assumption holds

Kesimpulan

  • Asumsi paralelisme penting untuk validitas model cumulative logit.

  • Menyederhanakan interpretasi karena efek prediktor konstan.

  • Jika tidak terpenuhi, gunakan model ordinal alternatif.

Log Linear Model

Analisis data kategorik penting dalam statistika terapan karena banyak fenomena nyata menghasilkan data dalam bentuk kategori, seperti jenis kelamin, pendidikan, atau diagnosis medis. Tiga pendekatan umum untuk menganalisis data ini adalah tabel kontingensi, model log-linier, dan regresi logistik, masing-masing dengan tujuan dan kelebihan berbeda.

Tabel kontingensi digunakan untuk eksplorasi awal hubungan antar variabel kategorik, misalnya untuk menghitung odds ratio atau statistik chi-kuadrat dalam studi efek obat terhadap penyakit.

Model log-linier cocok untuk menganalisis struktur asosiasi dalam tabel kontingensi yang kompleks. Model ini memodelkan frekuensi sel dengan distribusi Poisson dan memperlakukan semua variabel secara setara, tanpa membedakan variabel dependen dan independen.

Regresi logistik digunakan ketika ada satu variabel kategorik sebagai outcome dan beberapa variabel prediktor. Model ini memodelkan peluang suatu kejadian (logit) dan umum dipakai untuk prediksi, baik dalam studi observasional maupun eksperimental. Model ini juga memiliki versi untuk outcome dengan lebih dari dua kategori.

Kesimpulan: - Tabel Kontingensi → Deskriptif - Model Log-linier → Eksploratif - Regresi Logistik → Prediktif

Pemilihan metode tergantung pada tujuan analisis: apakah ingin mendeskripsikan, mengeksplorasi hubungan, atau memprediksi outcome. Kombinasi ketiganya sering digunakan untuk pemahaman yang menyeluruh terhadap data kategorik.

Tabel Kontingensi

Tabel kontingensi menyajikan jumlah frekuensi dari kombinasi kategori antar variabel.

Contoh :

table_data <- matrix(c(30, 20, 50, 70), nrow=2, 
       dimnames = list(Obat = c("Timolol", "Placebo"),
                       Serangan = c("Ya", "Tidak")))
table_data
##          Serangan
## Obat      Ya Tidak
##   Timolol 30    50
##   Placebo 20    70

Tabel kontingensi bersifat deskriptif dan tidak melibatkan pemodelan probabilitas.

Model Loglinear

Model loglinear memodelkan logaritma dari ekspektasi frekuensi sel dalam tabel kontingensi.

\[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \]

  • Cocok untuk menyelidiki asosiasi dan independensi antar variabel.
  • Tidak membedakan antara variabel respon dan penjelas.
  • Umumnya digunakan pada tabel multi-dimensi (2 arah, 3 arah, dst).
library(MASS) 
loglm(~ Obat * Serangan, data = table_data)
## Call:
## loglm(formula = ~Obat * Serangan, data = table_data)
## 
## Statistics:
##                  X^2 df P(> X^2)
## Likelihood Ratio   0  0        1
## Pearson            0  0        1

Model Regresi Logistik

Model regresi logistik biner:

\[ \log\left( \frac{p}{1 - p} \right) = \beta_0 + \beta_1 x \]

  • Digunakan jika ada variabel dependen kategorik (biasanya biner).

  • Bertujuan untuk memprediksi probabilitas suatu outcome.

  • Umumnya digunakan dalam studi observasional atau eksperimental.

Contoh :

data_glm <- data.frame(
  Serangan = c(1, 0, 1, 0),
  Obat = factor(c("Timolol", "Timolol", "Placebo", "Placebo")),
  Frek = c(30, 20, 50, 70)
)
model_logit <- glm(Serangan ~ Obat, weights = Frek, family = binomial, data = data_glm)
summary(model_logit)
## 
## Call:
## glm(formula = Serangan ~ Obat, family = binomial, data = data_glm, 
##     weights = Frek)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)  
## (Intercept)  -0.3365     0.1852  -1.817   0.0692 .
## ObatTimolol   0.7419     0.3430   2.163   0.0305 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 235.08  on 3  degrees of freedom
## Residual deviance: 230.31  on 2  degrees of freedom
## AIC: 234.31
## 
## Number of Fisher Scoring iterations: 4

Tabel Kontigensi dan Model Loglinier

Tabel kontingensi menyajikan frekuensi dari kombinasi kategori antar dua atau lebih variabel. Misal:

# Contoh tabel 2x2
matrix(c(30, 20, 50, 70), nrow=2,
       dimnames = list(Obat = c("Timolol", "Placebo"),
                       Serangan = c("Ya", "Tidak")))
##          Serangan
## Obat      Ya Tidak
##   Timolol 30    50
##   Placebo 20    70

Model log-linier untuk tabel I x J dapat dituliskan: \[ \log(\mu_{ij}) = \lambda + \lambda_i^{T} + \lambda_j^{R} + \lambda_{ij}^{TR} \]

Model Saturated

Model saturated atau model penuh menyertakan seluruh efek utama dan interaksi:

  • Cocok sempurna terhadap data
  • Tidak mengasumsikan independensi antar variabel Contoh formulasi untuk tabel 2x2:
# Data
library(MASS)
data <- matrix(c(35, 65, 45, 55), nrow=2, byrow=TRUE)
dimnames(data) <- list(Obat = c("Timolol", "Placebo"), Serangan = c("Ya", "Tidak"))
ftable(data)
##         Serangan Ya Tidak
## Obat                     
## Timolol          35    65
## Placebo          45    55

Model saturated dapat dipasang dengan loglm dari package {MASS}:

model_saturated <- loglm(~ Obat * Serangan, data = data)
summary(model_saturated)
## Formula:
## ~Obat * Serangan
## attr(,"variables")
## list(Obat, Serangan)
## attr(,"factors")
##          Obat Serangan Obat:Serangan
## Obat        1        0             1
## Serangan    0        1             1
## attr(,"term.labels")
## [1] "Obat"          "Serangan"      "Obat:Serangan"
## attr(,"order")
## [1] 1 1 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
## 
## Statistics:
##                  X^2 df P(> X^2)
## Likelihood Ratio   0  0        1
## Pearson            0  0        1

Model Independent

Model independen mengasumsikan bahwa tidak ada interaksi antara variabel: \[ \log(\mu_{ij}) = \lambda + \lambda_i^{T} + \lambda_j^{R} \]

Model ini menguji hipotesis bahwa variabel X dan Y saling independen.

model_indep <- loglm(~ Obat + Serangan, data = data)
summary(model_indep)
## Formula:
## ~Obat + Serangan
## attr(,"variables")
## list(Obat, Serangan)
## attr(,"factors")
##          Obat Serangan
## Obat        1        0
## Serangan    0        1
## attr(,"term.labels")
## [1] "Obat"     "Serangan"
## attr(,"order")
## [1] 1 1
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
## 
## Statistics:
##                       X^2 df  P(> X^2)
## Likelihood Ratio 2.087576  1 0.1485015
## Pearson          2.083333  1 0.1489147

Odds Ratio dan Interpretasi

Odds ratio untuk tabel 2x2:

\[ OR = \frac{n_{11} n_{22}}{n_{12} n_{21}} \]

Interpretasi nilai OR:

  • OR = 1: Tidak ada asosiasi
  • OR > 1: Asosiasi positif
  • OR < 1: Asosiasi negatif

Estimasi Parameter

Dalam model saturated:

  • Estimasi dilakukan dengan pembatasan seperti sum-to-zero
  • Estimasi parameter dilakukan dengan iterative proportional fitting (IPF)
# Estimasi odds ratio dan log-odds
logOR <- log((data[1,1] * data[2,2]) / (data[1,2] * data[2,1]))
logOR
## [1] -0.4183685

Model Lebih Sederhana dan Perbandingan Model

Perbandingan antar model dilakukan dengan menggunakan statistik deviance (G²) atau likelihood ratio test.

anova(model_indep, model_saturated)
## LR tests for hierarchical log-linear models
## 
## Model 1:
##  ~Obat + Serangan 
## Model 2:
##  ~Obat * Serangan 
## 
##           Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1   2.087576  1                                    
## Model 2   0.000000  0   2.087576         1         0.1485
## Saturated 0.000000  0   0.000000         0         1.0000

Studi Kasus: Kepercayaan terhadap Surga

Studi dari Agresti (2019) membahas hubungan antara kebahagiaan dan kepercayaan terhadap kehidupan akhirat.

data_survey <- matrix(c(32,190,
                        113,611,
                        51,326),
                      nrow = 3, byrow = TRUE,
                      dimnames = list(Kebahagiaan = c("Tidak", "Cukup", "Sangat"),
                                      Surga = c("Tidak Percaya", "Percaya")))
ftable(data_survey)
##             Surga Tidak Percaya Percaya
## Kebahagiaan                            
## Tidak                        32     190
## Cukup                       113     611
## Sangat                       51     326
loglm(~ Kebahagiaan + Surga, data = data_survey)
## Call:
## loglm(formula = ~Kebahagiaan + Surga, data = data_survey)
## 
## Statistics:
##                        X^2 df  P(> X^2)
## Likelihood Ratio 0.8911136  2 0.6404675
## Pearson          0.8836760  2 0.6428538

Rangkuman Model Log Linear 2 Arah

Model Log-Linear pada tabel kontingensi

Model log-linear adalah model yang digunakan untuk menganalisis hubungan antara dua atau lebih variabel kategorik yang disajikan dalam tabel kontingensi. Model ini mengasumsikan bahwa logaritma dari nilai ekspektasi frekuensi sel (μij ) dapat dinyatakan sebagai penjumlahan efek variabel dan (bila perlu) interaksinya. Untuk tabel 2x2: \[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \]

Perbedaan Utama Antara Model Log Linear dan Model Regresi Logistik

  • Model log-linear digunakan untuk memodelkan frekuensi (count) pada tabel kontingensi dan menguji asosiasi antar variabel kategorik, tanpa menganggap ada variabel respon dan prediktor.

  • Model regresi logistik digunakan untuk memodelkan probabilitas kejadian suatu outcome (biner) berdasarkan satu atau lebih prediktor (bisa kategorik maupun kontinu).

Estimasi Parameter Log Linear 2 Arah

Sistem Persamaan Model Log-Linear

\[\begin{align*} \log(\mu_{11}) &= \lambda + \lambda^A_1 + \lambda^B_1 + \lambda^{AB}_{11} \\ \log(\mu_{12}) &= \lambda + \lambda^A_1 + \lambda^B_2 + \lambda^{AB}_{12} \\ \log(\mu_{21}) &= \lambda + \lambda^A_2 + \lambda^B_1 + \lambda^{AB}_{21} \\ \log(\mu_{22}) &= \lambda + \lambda^A_2 + \lambda^B_2 + \lambda^{AB}_{22} \end{align*}\]

Constraint Sum-to-Zero

\[\begin{align*} \lambda_1^A + \lambda_2^A = 0 \\ \\ \lambda_1^B + \lambda_2^B = 0 \\ \\ \lambda_{11}^{AB} + \lambda_{12}^{AB} + \lambda_{21}^{AB} + \lambda_{22}^{AB} = 0 \end{align*}\]

Rumus Estimasi Parameter dengan Sum-to-Zero Constraint

\[ \lambda_1^A = \frac{1}{2} \left[ (\log \mu_{11} + \log \mu_{12}) - (\log \mu_{21} + \log \mu_{22}) \right] \] \[ \lambda_1^B = \frac{1}{2} \left[ (\log \mu_{11} + \log \mu_{21}) - (\log \mu_{12} + \log \mu_{22}) \right] \] \[ \lambda_{12}^{AB} = \frac{1}{4} \left[ \log \mu_{12} - \log \mu_{11} - \log \mu_{22} + \log \mu_{21} \right] \]

Analisis Data Tabel Kontingensi 2x2

Diberikan data:

Merokok Sakit Sehat
Ya 30 20
Tidak 10 40

Bentuk Model Log-Linear

Model log-linear pada tabel 2x2:

\[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \] dengan constraint sum-to-zero: \[ \sum_i \lambda_i^A = 0,\quad \sum_j \lambda_j^B = 0,\quad \sum_{i,j} \lambda_{ij}^{AB} = 0 \]

Estimasi Parameter Model (Manual, Sum-to-zero)

Misalkan:

  • A1 = Merokok (Ya), A2 = Tidak

  • B1 = Sakit, B2 = Sehat

Observasi:

  • \(n_{11}\) = 30, \(n_{12}\) = 20

  • \(n_{21}\) = 10, \(n_{22}\) = 40

\[\begin{align*} \log(\mu_{11}) &= \lambda + \lambda^A_1 + \lambda^B_1 + \lambda^{AB}_{11} \\ \\ \log(\mu_{12}) &= \lambda + \lambda^A_1 + \lambda^B_2 + \lambda^{AB}_{12} \\ \\ \log(\mu_{21}) &= \lambda + \lambda^A_2 + \lambda^B_1 + \lambda^{AB}_{21} \\ \\ \log(\mu_{22}) &= \lambda + \lambda^A_2 + \lambda^B_2 + \lambda^{AB}_{22} \\ \end{align*} \] Constraint sum-to-zero:

\[\begin{align*} \lambda_1^A + \lambda_2^A = 0 \\ \\ \lambda_1^B + \lambda_2^B = 0 \\ \\ \lambda_{11}^{AB} + \lambda_{12}^{AB} + \lambda_{21}^{AB} + \lambda_{22}^{AB} = 0 \end{align*}\]

Langkah-langkah:

  1. Hitung rata-rata log frekuensi sel:

\[\begin{align*} \lambda &= \frac{1}{4} \sum_{i=1}^{2} \sum_{j=1}^{2} \log(n_{ij}) \\ &= \frac{1}{4} [\log(30) + \log(20) + \log(10) + \log(40)] \\ &= 3.0971 \end{align*} \]

  1. Efek Utama A (Merokok):

\[\begin{align*} \lambda_1^A &= \frac{1}{2} [(\log(30) + \log(20)) - (\log(10) + \log(40))] \\ \\ &= \frac{1}{2} [(3.4012 + 2.9957) - (2.3026 + 3.6889)] \\ \\ &= \frac{1}{2} (6.3969 - 5.9915) \\ \\ &= \frac{1}{2} (0.4054) \\ \\ &= 0.2027 \\ \\ \lambda_2^A &= -0.2027 \\ \end{align*}\]

  1. Efek Utama B (Status):

\[\begin{align*} \\ \lambda_1^B &= \frac{1}{2} [(\log(30) + \log(10)) - (\log(20) + \log(40))] \\ \\ &= \frac{1}{2} [(3.4012 + 2.3026) - (2.9957 + 3.6889)] \\ \\ &= \frac{1}{2} (5.7038 - 6.6846) \\ \\ &= \frac{1}{2} (-0.9808) \\ \\ &= -0.4904 \\ \\ \lambda_2^B &= +0.4904 \\ \end{align*}\]

  1. Efek Interaksi:

\[\begin{align*} \lambda_{11}^{AB} &= \frac{1}{4} [(\log(30) - \log(20)) - \log(10) + \log(40)] \\ \\ &= \frac{1}{4} [(3.4012 - 2.9957 - 2.3026 + 3.6889)] \\ \\ &= \frac{1}{4} (-1.8971 + 3.6889) \\ \\ &= \frac{1}{4} (1.7918) \\ \\ &= 0.4479 \\ \\ \lambda_{12}^{AB} &= -0.4479 \\ \\ \lambda_{21}^{AB} &= -0.4479 \\ \\ \lambda_{22}^{AB} &= 0.4479 \\ \end{align*}\]

Ringkasan Parameter

  • \(\lambda = 3.0971\)
  • \(\lambda_1^A = 0.2027, \lambda_2^A = -0.2027\)
  • \(\lambda_1^B = -0.4904, \lambda_2^B = 0.4904\)
  • \(\lambda_{11}^{AB} = 0.4479, \lambda_{12}^{AB} = -0.4479, \lambda_{21}^{AB} = -0.4479, \lambda_{22}^{AB} = 0.4479\)

Hitung Odds Ratio dan Interval Kepercayaan

\[ OR = \frac{n_{11} n_{22}}{n_{12} n_{21}} = \frac{30 \times 40}{20 \times 10} = \frac{1200}{200} = 6 \] Log odds ratio: \[ \log(OR) = \log(6) = 1.7918 \] Standard error (SE):

\[ SE = \sqrt{ \frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}} } \\ \\ = \sqrt{ \frac{1}{30} + \frac{1}{20} + \frac{1}{10} + \frac{1}{40} }\\ \\ = \sqrt{0.0333 + 0.05 + 0.1 + 0.025}\\ \\ = \sqrt{0.2083}\\ \\ = 0.4564 \\ \]

95% Confidence Interval for log(OR):

\[\begin{align*} \\ \log(OR) \pm 1.96 \times SE &= 1.7918 \pm 1.96 \times 0.4564 \\ \\ &= 1.7918 \pm 0.895 \\ \\ &= (1.7918 - 0.895,\ 1.7918 + 0.895) \\ \\ &= (0.8968,\ 2.6868) \\ \end{align*}\]

Back-transform to get CI for OR:

\[\begin{align*} \\ Lower &= \exp(0.8968) = 2.452\\ \\ Upper &= \exp(2.6868) = 14.68 \\ \end{align*}\]

Jadi, OR=6 dengan taraf signifikansi 95% memiliki interval kepercayaan dalam rentang 2.452 - 14.68

Fitting Model Log-Linear dengan R

# Data 2x2
tabel <- matrix(c(30, 20, 10, 40), nrow = 2, byrow = TRUE)
colnames(tabel) <- c("Sakit", "Sehat")
rownames(tabel) <- c("Ya", "Tidak")
tabel
##       Sakit Sehat
## Ya       30    20
## Tidak    10    40
#Ubah menjadi bentuk tabel untuk glm
data <- as.data.frame(as.table(tabel))
colnames(data) <- c("Merokok", "Status", "Freq")
data
# Model tanpa interaksi
fit_no_inter <- glm(Freq ~ Merokok + Status, family = poisson, data = data)
summary(fit_no_inter)
## 
## Call:
## glm(formula = Freq ~ Merokok + Status, family = poisson, data = data)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  2.996e+00  1.871e-01  16.013   <2e-16 ***
## MerokokTidak 3.892e-10  2.000e-01   0.000    1.000    
## StatusSehat  4.055e-01  2.041e-01   1.986    0.047 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 21.288  on 3  degrees of freedom
## Residual deviance: 17.261  on 1  degrees of freedom
## AIC: 43.036
## 
## Number of Fisher Scoring iterations: 4
# Model dengan interaksi
fit_inter <- glm(Freq ~ Merokok * Status, family = poisson, data = data)
summary(fit_inter)
## 
## Call:
## glm(formula = Freq ~ Merokok * Status, family = poisson, data = data)
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                3.4012     0.1826  18.629  < 2e-16 ***
## MerokokTidak              -1.0986     0.3651  -3.009  0.00262 ** 
## StatusSehat               -0.4055     0.2887  -1.405  0.16015    
## MerokokTidak:StatusSehat   1.7918     0.4564   3.926 8.65e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2.1288e+01  on 3  degrees of freedom
## Residual deviance: 3.9968e-15  on 0  degrees of freedom
## AIC: 27.775
## 
## Number of Fisher Scoring iterations: 3

Interpretasi Parameter

  • Parameter Utama (Intercept) menunjukkan rata-rata log frekuensi sel
  • Efek “Merokok” dan “Status”menunjukkan perbedaan log frekuensi antar kategori
  • Interaksi signifika menunjukkan adanya asosiasi antara efek Merokok dan Status

Nilai \(\log(6)=1.7918\) sama dengan efek interaksi pada output R

Analisis Data Tabel Kontingensi 2x3

Suatu survei dilakukan untuk mengetahui hubungan antara Jenis Kelamin (Laki-laki/Perempuan) dan Kategori BMI (Kurus/Normal/Gemuk):

Kurus Normal Gemuk
Laki-laki 12 20 8
Perempuan 18 24 10

Bentuk Model Log-Linear untuk Tabel 2x3

Bentuk umum model log-linear untuk tabel 2x3 (dengan sum-to-zero constraint): \[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \] dengan: - \(\mu_{ij}\): ekspektasi frekuensi pada baris ke-i, kolom ke-j - A: Jenis Kelamin (i=1:Laki-laki, i=2:Perempuan) - B: Kategori BMI (j=1:Kurus, j=2:Normal, j=3:Gemuk) - Constraint: \[ \sum_i \lambda_i^A = 0,\quad \sum_j \lambda_j^B = 0,\quad \sum_{i} \lambda_{ij}^{AB} = 0, dan \sum_{j}\lambda_{ij}^{AB}=0 \] Secara eksplisit:

\[\begin{align*} \\ \log(\mu_{ij}) &= \lambda \\ \\ &+ \lambda^A_1(Laki-laki), \lambda_2^A(Perempuan)\\ &+ \lambda^B_1(Kurus), \lambda^B_2(Normal), \lambda^B_3(Gemuk)\\ \\ &+ \lambda^{AB}_{ij}(interaksi jika ada) \\ \end{align*}\]

Fitting Model Log-Linear di R

# Membuat data frame dari tabel
tabel2x3 <- matrix(c(12, 20, 8, 18, 24, 10), nrow = 2, byrow = TRUE)
colnames(tabel2x3) <- c("Kurus", "Normal", "Gemuk")
rownames(tabel2x3) <- c("Laki-laki", "Perempuan")
tabel2x3
##           Kurus Normal Gemuk
## Laki-laki    12     20     8
## Perempuan    18     24    10
# Ubah menjadi data.frame untuk glm
data2x3 <- as.data.frame(as.table(tabel2x3))
colnames(data2x3) <- c("JenisKelamin", "BMI", "Freq")
data2x3
# Model log-linear tanpa interaksi (asumsi independen)
fit_no_inter <- glm(Freq ~ JenisKelamin + BMI, family = poisson, data = data2x3)
summary(fit_no_inter)
## 
## Call:
## glm(formula = Freq ~ JenisKelamin + BMI, family = poisson, data = data2x3)
## 
## Coefficients:
##                       Estimate Std. Error z value Pr(>|z|)    
## (Intercept)             2.5683     0.2179  11.789   <2e-16 ***
## JenisKelaminPerempuan   0.2624     0.2103   1.248   0.2122    
## BMINormal               0.3830     0.2368   1.618   0.1058    
## BMIGemuk               -0.5108     0.2981  -1.713   0.0866 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 13.06443  on 5  degrees of freedom
## Residual deviance:  0.22527  on 2  degrees of freedom
## AIC: 35.26
## 
## Number of Fisher Scoring iterations: 3
# Model log-linear dengan interaksi (untuk cek asosiasi)
fit_inter <- glm(Freq ~ JenisKelamin * BMI, family = poisson, data = data2x3)
summary(fit_inter)
## 
## Call:
## glm(formula = Freq ~ JenisKelamin * BMI, family = poisson, data = data2x3)
## 
## Coefficients:
##                                 Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       2.4849     0.2887   8.608   <2e-16 ***
## JenisKelaminPerempuan             0.4055     0.3727   1.088    0.277    
## BMINormal                         0.5108     0.3651   1.399    0.162    
## BMIGemuk                         -0.4055     0.4564  -0.888    0.374    
## JenisKelaminPerempuan:BMINormal  -0.2231     0.4802  -0.465    0.642    
## JenisKelaminPerempuan:BMIGemuk   -0.1823     0.6032  -0.302    0.762    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance:  1.3064e+01  on 5  degrees of freedom
## Residual deviance: -9.0719e-30  on 0  degrees of freedom
## AIC: 39.034
## 
## Number of Fisher Scoring iterations: 3

Interpretasi:

Model tanpa interaksi:

  • Jika deviance tidak signifikan, maka Jenis Kelamin dan BMI independen.

  • Intercept: log frekuensi pada kategori referensi (Laki-laki, Kurus)

  • Koefisien JenisKelaminPerempuan: perbedaan log-frekuensi antara Perempuan vs Laki-laki (pada Kurus)

  • Koefisien BMI: perbedaan log-frekuensi kategori BMI (Normal/Gemuk) terhadap Kurus (pada Laki-laki)

Model dengan interaksi:

  • Jika koefisien interaksi signifikan, berarti ada hubungan/asosiasi antara Jenis Kelamin dan BMI. Artinya distribusi BMI berbeda antara Laki-laki dan Perempuan.

Contoh interpretasi hasil (misal)

  • Jika koefisien JenisKelaminPerempuan negatif: proporsi Perempuan pada kategori referensi lebih kecil dibanding Laki-laki.

  • Jika koefisien BMI_Normal positif: kemungkinan seseorang Normal lebih tinggi daripada Kurus (pada Laki-laki).

  • Jika model interaksi signifikan, pola distribusi BMI berbeda antara Laki-laki dan Perempuan.

Model Log Linear Tiga Arah

Pada pembahasan sebelumnya, kita telah memahami bahwa salah satu tujuan utama dari penyusunan model log linear adalah untuk mengestimasi parameter-parameter yang menjelaskan hubungan di antara variabel-variabel kategorik.

Pada materi kali ini, kita akan membahas model log linear yang lebih kompleks, yaitu model log linear untuk tabel kontingensi tiga arah. Model ini melibatkan tiga variabel kategorik, sehingga kemungkinan interaksi yang dapat terjadi di dalam model pun menjadi lebih banyak. Dalam konteks ini, interaksi paling tinggi yang dapat dimodelkan adalah interaksi tiga arah, yaitu interaksi yang melibatkan ketiga variabel secara bersamaan.

Model Log-Linear untuk Tabel Tiga Arah

Model log-linear yang melibatkan tiga variabel kategorik (misal: X, Y, dan Z) dapat dibangun dalam berbagai bentuk model, tergantung pada tingkat interaksi yang ingin dimasukkan. Berikut adalah beberapa alternatif model log-linear yang umum digunakan: 1. Model Saturated \[ \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 memuat semua kemungkinan interaksi, termasuk interaksi tiga arah (X, Y, dan Z).

  1. Model Homogen \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \] Model ini hanya mengakomodasi interaksi dua arah antar variabel tanpa memasukkan interaksi tiga arah.

  2. Model Conditional

    • Conditional pada X: \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \] Memuat interaksi X dengan Y dan X dengan Z.

    • Conditional pada Y: \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \] Memuat interaksi Y dengan X dan Y dengan Z.

    • Conditional pada Z: \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \] Memuat interaksi Z dengan X dan Z dengan Y.

  3. Model Joint Independence

    • 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} \]

  4. Model Tanpa Interaksi \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k \] Model ini hanya memasukkan efek utama tanpa interaksi antar variabel.

Pengujian Interaksi dalam Model Log-Linear 3 Arah

Dalam analisis model log-linear tiga arah, pengujian interaksi dilakukan untuk mengetahui ada atau tidaknya interaksi antar variabel. Pengujian ini dilakukan secara bertahap, dimulai dari tingkat interaksi tertinggi ke yang lebih rendah. Untuk model log-linear dengan tiga peubah (X, Y, dan Z), tahapan pengujian meliputi:

  1. Pengujian Interaksi Tiga Arah (XYZ):
    • Bandingkan model saturated dengan model homogenous.
  2. Pengujian Interaksi Dua Arah (XY, XZ, YZ):

    • Bandingkan model homogenous dengan model conditional.

    • Bandingkan model conditional dengan model joint independence.

    • Bandingkan model joint independence dengan model tanpa interaksi.

Setiap tahapan pengujian dilakukan untuk menilai kecocokan model dan menentukan struktur interaksi mana yang paling sesuai dengan data yang diamati.

Soal Praktikum

Tabel berikut menyajikan data dari survei General Social Survey (GSS) tahun 1994 mengenai jenis kelamin responden, tingkat fundamentalisme, dan sikap terhadap hukuman mati untuk kasus pembunuhan. Susun dan interpretasikan model log-linear paling sederhana (paling parsimonious) untuk data ini. Jelaskan proses yang Anda lakukan dalam menentukan model terbaik serta asosiasi apa saja yang teridentifikasi. Tunjukkan juga bagaimana nilai yang diprediksi dari model menggambarkan asosiasi tersebut.

Fundamentalism Jenis Kelamin Mendukung Menolak Total
Fundamentalist Laki-laki 128 32 160
Fundamentalist Perempuan 123 73 196
Fundamentalist Total 251 105 356
Moderate Laki-laki 182 56 238
Moderate Perempuan 168 105 273
Moderate Total 350 161 511
Liberal Laki-laki 119 49 168
Liberal Perempuan 111 70 181
Liberal Total 230 119 349

Keterangan:

  • Fundamentalisme: Fundamentalist, Moderate, Liberal

  • Jenis Kelamin: Laki-laki, Perempuan

  • Sikap: Mendukung (Favor), Menolak (Oppose) hukuman mati

Analisis Log-Linear untuk Tabel Tiga Arah

library("epitools")
library("DescTools")
## Warning: package 'DescTools' was built under R version 4.3.3
library("lawstat")
## Warning: package 'lawstat' was built under R version 4.3.3
# Input data sesuai tabel praktikum
z.fund <- factor(rep(c("1fund", "2mod", "3lib"), each = 4))
x.sex  <- factor(rep(c("1M", "2F"), each = 2, times = 3))
y.fav  <- factor(rep(c("1fav", "2opp"), times = 6))
counts <- c(128, 32, 123, 73, 182, 56, 168, 105, 119, 49, 111, 70)

data <- data.frame(
  Fundamentalisme = z.fund,
  Jenis_Kelamin   = x.sex,
  Sikap           = y.fav,
  Frekuensi       = counts

)
data

Membentuk Tabel Kontingensi 3 Arah

table3d <- xtabs(Frekuensi ~ Fundamentalisme + Jenis_Kelamin + Sikap, data = data)
ftable(table3d)
##                               Sikap 1fav 2opp
## Fundamentalisme Jenis_Kelamin                
## 1fund           1M                   128   32
##                 2F                   123   73
## 2mod            1M                   182   56
##                 2F                   168  105
## 3lib            1M                   119   49
##                 2F                   111   70

Analisis Log-Linear: Tahap Pemodelan Kita akan memodelkan tabel ini menggunakan beberapa model log-linear dan membandingkan kecocokan model (parsimonious model):

Uji Model Interaksi Tiga Arah (Saturated VS Homogenous)

Penentuan Kategori Referensi

##=============================##
# Penentuan kategori reference
##=============================##
x.sex  <- relevel(x.sex, ref = "2F")
y.fav  <- relevel(y.fav, ref = "2opp")
z.fund <- relevel(z.fund, ref = "3lib")

Model Saturated memasukkan semua interaksi hingga 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 saturated
model_saturated <- glm(counts ~ x.sex + y.fav + z.fund +
             x.sex*y.fav + x.sex*z.fund + y.fav*z.fund +
             x.sex*y.fav*z.fund,
             family = poisson(link = "log"))
summary(model_saturated)
## 
## Call:
## glm(formula = counts ~ x.sex + y.fav + z.fund + x.sex * y.fav + 
##     x.sex * z.fund + y.fav * z.fund + x.sex * y.fav * z.fund, 
##     family = poisson(link = "log"))
## 
## Coefficients:
##                                Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                    4.248495   0.119523  35.545  < 2e-16 ***
## x.sex1M                       -0.356675   0.186263  -1.915  0.05551 .  
## y.fav1fav                      0.461035   0.152626   3.021  0.00252 ** 
## z.fund1fund                    0.041964   0.167285   0.251  0.80193    
## z.fund2mod                     0.405465   0.154303   2.628  0.00860 ** 
## x.sex1M:y.fav1fav              0.426268   0.228268   1.867  0.06185 .  
## x.sex1M:z.fund1fund           -0.468049   0.282210  -1.659  0.09721 .  
## x.sex1M:z.fund2mod            -0.271934   0.249148  -1.091  0.27507    
## y.fav1fav:z.fund1fund          0.060690   0.212423   0.286  0.77511    
## y.fav1fav:z.fund2mod           0.008969   0.196903   0.046  0.96367    
## x.sex1M:y.fav1fav:z.fund1fund  0.438301   0.336151   1.304  0.19227    
## x.sex1M:y.fav1fav:z.fund2mod   0.282383   0.301553   0.936  0.34905    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2.4536e+02  on 11  degrees of freedom
## Residual deviance: 5.9952e-15  on  0  degrees of freedom
## AIC: 100.14
## 
## Number of Fisher Scoring iterations: 3
exp(model_saturated$coefficients)
##                   (Intercept)                       x.sex1M 
##                    70.0000000                     0.7000000 
##                     y.fav1fav                   z.fund1fund 
##                     1.5857143                     1.0428571 
##                    z.fund2mod             x.sex1M:y.fav1fav 
##                     1.5000000                     1.5315315 
##           x.sex1M:z.fund1fund            x.sex1M:z.fund2mod 
##                     0.6262231                     0.7619048 
##         y.fav1fav:z.fund1fund          y.fav1fav:z.fund2mod 
##                     1.0625694                     1.0090090 
## x.sex1M:y.fav1fav:z.fund1fund  x.sex1M:y.fav1fav:z.fund2mod 
##                     1.5500717                     1.3262868

Model yang digunakan adalah model log-linear saturated dengan semua efek utama, interaksi dua arah, dan interaksi tiga arah. Model ini memodelkan hubungan antara jenis kelamin (x.sex), sikap terhadap hukuman mati (y.fav), dan tingkat fundamentalisme (z.fund) terhadap frekuensi responden. Dengan hasil estimasi koefisien nya sebagai berikut:

Parameter Estimate Std. Error z value Pr(>z)
(Intercept) 4.25 0.12 35.55 <2e-16***
x.sex1M -0.36 0.19 -1.92 0.055 .
y.fav1fav 0.46 0.15 3.02 0.0025 **
z.fund1fund 0.04 0.17 0.25 0.80
z.fund2mod 0.41 0.15 2.63 0.0086 **
x.sex1M:y.fav1fav 0.43 0.23 1.87 0.062 .
x.sex1M:z.fund1fund -0.47 0.28 -1.66 0.097 .
x.sex1M:z.fund2mod -0.27 0.25 -1.09 0.28
y.fav1fav:z.fund1fund 0.06 0.21 0.29 0.78
y.fav1fav:z.fund2mod 0.01 0.20 0.05 0.96
x.sex1M:y.fav1fav:z.fund1fund 0.44 0.34 1.30 0.19
x.sex1M:y.fav1fav:z.fund2mod 0.28 0.30 0.94 0.35

Interpretasi Koefisien:

Goodness-of-Fit:

Kesimpulan:

Catatan Interpretasi:

Model Homogenous

Model log-linear homogenous memasukkan semua efek utama dan semua interaksi dua arah, tanpa interaksi tiga arah. Secara matematis, model ini dapat dituliskan sebagai berikut: \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]

# Homogenous Model
model_homogenous <- glm(counts ~ x.sex + y.fav + z.fund +
              x.sex*y.fav + x.sex*z.fund + y.fav*z.fund,
              family = poisson(link = "log"))
summary(model_homogenous)
## 
## Call:
## glm(formula = counts ~ x.sex + y.fav + z.fund + x.sex * y.fav + 
##     x.sex * z.fund + y.fav * z.fund, family = poisson(link = "log"))
## 
## Coefficients:
##                       Estimate Std. Error z value Pr(>|z|)    
## (Intercept)            4.31096    0.10522  40.972  < 2e-16 ***
## x.sex1M               -0.51575    0.13814  -3.733 0.000189 ***
## y.fav1fav              0.35707    0.12658   2.821 0.004788 ** 
## z.fund1fund           -0.06762    0.14452  -0.468 0.639854    
## z.fund2mod             0.33196    0.13142   2.526 0.011540 *  
## x.sex1M:y.fav1fav      0.66406    0.12728   5.217 1.81e-07 ***
## x.sex1M:z.fund1fund   -0.16201    0.15300  -1.059 0.289649    
## x.sex1M:z.fund2mod    -0.08146    0.14079  -0.579 0.562887    
## y.fav1fav:z.fund1fund  0.23873    0.16402   1.455 0.145551    
## y.fav1fav:z.fund2mod   0.13081    0.14951   0.875 0.381614    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 245.361  on 11  degrees of freedom
## Residual deviance:   1.798  on  2  degrees of freedom
## AIC: 97.934
## 
## Number of Fisher Scoring iterations: 3

Uji Hipotesis: Apakah Ada Interaksi Tiga Arah? (Saturated vs Homogenous)

Pengujian ini menggunakan residual deviance dari kedua model (saturated dan homogenous).

Langkah-Langkah Pengujian:

  1. Hipotesis:

  • H0: Tidak ada interaksi tiga arah (model homogenous sudah cukup)

  • H1: Ada interaksi tiga arah (model saturated diperlukan)

  1. Hitung Selisih Deviance
# Deviance antar model
Deviance.model <- model_homogenous$deviance - model_saturated$deviance
Deviance.model
## [1] 1.797977
  1. Hitung Derajat Bebas
# Derajat bebas = db model homogenous - db model saturated
derajat.bebas <- (model_homogenous$df.residual - model_saturated$df.residual)
derajat.bebas
## [1] 2
  1. Chi-Square Tabel (\(\alpha\)=0.05)
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel
## [1] 5.991465
  1. Keputusan Uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan
## [1] "Terima H0"

Interpretasi Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara jenis kelamin, fundamentalisme, dan pendapat mengenai hukuman mati.

Catatan:

  • Model pengurang adalah model yang lebih lengkap (lebih banyak parameter, df lebih kecil), yaitu model saturated

  • Derajat bebas dihitung dari selisih derajat bebas model homogenous dan saturated.

  • Keputusan berdasarkan perbandingan deviance model dengan chi-square tabel.

Rangkuman

Pengujian ada tidaknya interaksi tiga arah (Saturated Model vs Homogenous Model)

  • Hipotesis

    • H0: \(\lambda_{ijk}^{XYZ}=0\) (Tidak ada interaksi tiga arah; model yang terbentuk adalah model homogenous) -

    • H0: \(\lambda_{ijk}^{XYZ}\neq0\) (Tidak ada interaksi tiga arah; model yang terbentuk adalah model homogenous)

  • Tingkat Signifikansi

    \(\alpha=5\%\)

  • Statistik Uji

    • \(\Delta\)Deviance = Deviance model homogenous - Deviance model saturated = 1.798 - 0.00 = 1.798

    • \(db\) = \(db\) model homogenous - \(db\) model saturated = 2 - 0 = 2

  • Daerah Penolakan

    Tolak H0 jika \(\Delta\)Deviance \(\>>\) \(\chi_{0.05,db}^2 = \chi_{0.05,2}^2 = 5.991\)

  • Keputusan

    Karena 1.798 < 5.991, maka terima H0

  • Interpretasi

    Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara jenis kelamin, fundamentalisme, dan pendapat mengenai hukuman mati.

Uji Model Interaksi Dua Arah (Homogenous vs Conditional on X)

Model log-linear conditional pada X memasukkan efek utama dan interaksi dua arah antara X dengan Y dan X dengan Z, tanpa interaksi antara Y dengan Z maupun interaksi tiga arah. \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \]

# Conditional Association on X
model_conditional_X <- glm(counts ~ x.sex + y.fav + z.fund +
              x.sex*y.fav + x.sex*z.fund,
              family = poisson(link = "log"))
summary(model_conditional_X)
## 
## Call:
## glm(formula = counts ~ x.sex + y.fav + z.fund + x.sex * y.fav + 
##     x.sex * z.fund, family = poisson(link = "log"))
## 
## Coefficients:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          4.23495    0.08955  47.293  < 2e-16 ***
## x.sex1M             -0.52960    0.13966  -3.792 0.000149 ***
## y.fav1fav            0.48302    0.08075   5.982 2.20e-09 ***
## z.fund1fund          0.07962    0.10309   0.772 0.439916    
## z.fund2mod           0.41097    0.09585   4.288 1.81e-05 ***
## x.sex1M:y.fav1fav    0.65845    0.12708   5.181 2.20e-07 ***
## x.sex1M:z.fund1fund -0.12841    0.15109  -0.850 0.395405    
## x.sex1M:z.fund2mod  -0.06267    0.13908  -0.451 0.652274    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 245.3612  on 11  degrees of freedom
## Residual deviance:   3.9303  on  4  degrees of freedom
## AIC: 96.067
## 
## Number of Fisher Scoring iterations: 4

Pengujian Ada Tidaknya Interaksi Antara Y dan Z (Homogenous Model vs Conditional Association on X)

  • Hipotesis

    • H0: \(\lambda_{jk}^{YZ}=0\) (Tidak ada interaksi antara pendapat hukuman mati (Y) dan fundamentalisme (Z)) -

    • H0: \(\lambda_{jk}^{YZ}\neq0\) (Ada interaksi antara pendapat hukuman mati (Y) dan fundamentalisme (Z))

  • Tingkat Signifikansi

    \(\alpha=5\%\)

  • Statistik Uji

    • \(\Delta\)Deviance = Deviance model conditional on X - Deviance model homogenous = 3.903 - 1.798 = 2.132

    • \(db\) = \(db\) model homogenous - \(db\) model saturated = 4 - 2 = 2

  • Daerah Penolakan

    Tolak H0 jika \(\Delta\)Deviance \(\>>\) \(\chi_{0.05,db}^2 = \chi_{0.05,2}^2 = 5.991\)

  • Keputusan

    Karena 2.132 < 5.991, maka terima H0

  • Interpretasi

    Dengan taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi antara pendapat tentang hukuman mati dan fundamentalisme. Dengan kata lain, model yang terbentuk adalah model tanpa parameter \(\lambda_{jk}^{YZ}\)

  1. Hitung Selisih Deviance
# Pengujian hipotesis

# Deviance of Model
Deviance.model <- model_conditional_X$deviance - model_homogenous$deviance   # model_conditional_X: conditional on X, model_homogenous: homogenous
Deviance.model
## [1] 2.132302
  1. Hitung Derajat Bebas
# Chi Square tabel dengan alpha = 0.05
derajat.bebas <- (4 - 2)
derajat.bebas
## [1] 2
  1. Chi-Square Tabel
chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 5.991465
  1. Keputusan Uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Terima"

Interpretasi Karena nilai Deviance.model = 2.13 lebih kecil dari nilai kritis chi-square tabel = 5.99 (dengan df = 2, alpha = 0.05), maka keputusan uji adalah “Terima”.

Pada taraf nyata 5%, belum cukup bukti untuk menolak H0, atau dengan kata lain tidak ada interaksi antara pendapat mengenai hukuman mati (Y) dan fundamentalisme (Z). Model yang terbentuk cukup hanya sampai dua interaksi dengan X (conditional on X), sehingga interaksi Y*Z tidak signifikan secara statistik.

Uji Model Interaksi Dua Arah (Homogenous vs Conditional on Y)

Model log-linear conditional pada Y memasukkan efek utama dan interaksi dua arah antara X dengan Y dan Y dengan Z, tanpa interaksi antara X dengan Z maupun interaksi tiga arah. \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \]

# Conditional Association on Y
model_conditional_Y <- glm(counts ~ x.sex + y.fav + z.fund +
              x.sex*y.fav + y.fav*z.fund,
              family = poisson(link = "log"))
summary(model_conditional_Y)
## 
## Call:
## glm(formula = counts ~ x.sex + y.fav + z.fund + x.sex * y.fav + 
##     y.fav * z.fund, family = poisson(link = "log"))
## 
## Coefficients:
##                       Estimate Std. Error z value Pr(>|z|)    
## (Intercept)            4.33931    0.09919  43.748  < 2e-16 ***
## x.sex1M               -0.59345    0.10645  -5.575 2.48e-08 ***
## y.fav1fav              0.37259    0.12438   2.996  0.00274 ** 
## z.fund1fund           -0.12516    0.13389  -0.935  0.34989    
## z.fund2mod             0.30228    0.12089   2.500  0.01240 *  
## x.sex1M:y.fav1fav      0.65845    0.12708   5.181 2.20e-07 ***
## y.fav1fav:z.fund1fund  0.21254    0.16205   1.312  0.18966    
## y.fav1fav:z.fund2mod   0.11757    0.14771   0.796  0.42606    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 245.3612  on 11  degrees of freedom
## Residual deviance:   2.9203  on  4  degrees of freedom
## AIC: 95.057
## 
## Number of Fisher Scoring iterations: 4

Pengujian Ada Tidaknya Interaksi Antara X dan Z (Homogenous Model vs Conditional Association on Y)

  • Hipotesis

    • H0: \(\lambda_{ik}^{XZ}=0\) (Tidak ada interaksi antara jenis kelamin (X) dan fundamentalisme (Z)) -

    • H0: \(\lambda_{ik}^{XZ}\neq0\) (Ada interaksi antara jenis kelamin

      1. dan fundamentalisme (Z))

  • Tingkat Signifikansi

    \(\alpha=5\%\)

  • Statistik Uji

    • \(\Delta\)Deviance = Deviance model conditional on Y - Deviance model homogenous = 2.9203 - 1.798 = 1.1223

    • \(db\) = \(db\) model homogenous - \(db\) model saturated = 4 - 2 = 2

  • Daerah Penolakan

    Tolak H0 jika \(\Delta\)Deviance \(\>>\) \(\chi_{0.05,db}^2 = \chi_{0.05,2}^2 = 5.991\)

  • Keputusan

    Karena 1.1223 < 5.991, maka terima H0

  • Interpretasi

    Dengan taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi antara jenis kelamin dan fundamentalisme. Dengan kata lain, model yang terbentuk adalah model tanpa parameter \(\lambda_{ik}^{XZ}\)

  1. Hitung Selisih Deviance
# Deviance of Model
Deviance.model <- model_conditional_Y$deviance - model_homogenous$deviance   # model_conditional_Y: conditional on Y, model_homogenous: homogenous
Deviance.model
## [1] 1.122315
  1. Hitung Derajat Bebas
# Chi Square tabel dengan alpha = 0.05
derajat.bebas <- (4 - 2)
derajat.bebas
## [1] 2
  1. Chi-Square Tabel
chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 5.991465
  1. Keputusan Uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Terima"

Interpretasi Karena nilai Deviance.model = 1.12 lebih kecil dari nilai kritis chi-square tabel = 5.99 (dengan df = 2, alpha = 0.05), maka keputusan uji adalah “Terima”.

Pada taraf nyata 5%, belum cukup bukti untuk menolak H0, atau dengan kata lain tidak ada interaksi antara jenis kelamin (X) dan fundamentalisme (Z). Model yang terbentuk cukup hanya sampai dua interaksi dengan Y (conditional on Y), sehingga interaksi X*Z tidak signifikan secara statistik.

Uji Model Interaksi Dua Arah (Homogenous vs Conditional on Z)

Model log-linear conditional pada Z memasukkan efek utama dan interaksi dua arah antara X dengan Z dan Y dengan Z, tanpa interaksi antara X dengan Y maupun interaksi tiga arah. \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]

# Conditional Association on Z
model_conditional_Z <- glm(counts ~ x.sex + y.fav + z.fund +
              x.sex*z.fund + y.fav*z.fund,
              family = poisson(link = "log"))
summary(model_conditional_Z)
## 
## Call:
## glm(formula = counts ~ x.sex + y.fav + z.fund + x.sex * z.fund + 
##     y.fav * z.fund, family = poisson(link = "log"))
## 
## Coefficients:
##                       Estimate Std. Error z value Pr(>|z|)    
## (Intercept)            4.12255    0.10518  39.195  < 2e-16 ***
## x.sex1M               -0.07453    0.10713  -0.696    0.487    
## y.fav1fav              0.65896    0.11292   5.836 5.36e-09 ***
## z.fund1fund           -0.06540    0.15126  -0.432    0.665    
## z.fund2mod             0.33196    0.13777   2.410    0.016 *  
## x.sex1M:z.fund1fund   -0.12841    0.15109  -0.850    0.395    
## x.sex1M:z.fund2mod    -0.06267    0.13908  -0.451    0.652    
## y.fav1fav:z.fund1fund  0.21254    0.16205   1.312    0.190    
## y.fav1fav:z.fund2mod   0.11757    0.14771   0.796    0.426    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 245.361  on 11  degrees of freedom
## Residual deviance:  29.729  on  3  degrees of freedom
## AIC: 123.87
## 
## Number of Fisher Scoring iterations: 4

Pengujian Ada Tidaknya Interaksi Antara X dan Y (Homogenous Model vs Conditional Association on Z)

  • Hipotesis

    • H0: \(\lambda_{ij}^{XY}=0\) (Tidak ada interaksi antara jenis kelamin (X) dan pendapat tentang hukuman mati (Y)) -

    • H0: \(\lambda_{ij}^{XY}\neq0\) (Ada interaksi antara jenis kelamin

      1. dan pendapat tentang hukuman mati (Y))

  • Tingkat Signifikansi

    \(\alpha=5\%\)

  • Statistik Uji

    • \(\Delta\)Deviance = Deviance model conditional on Y - Deviance model homogenous = 29.729 - 1.798 = 27.931

    • \(db\) = \(db\) model homogenous - \(db\) model saturated = 3 - 2 = 1

  • Daerah Penolakan

    Tolak H0 jika \(\Delta\)Deviance \(\>>\) \(\chi_{0.05,db}^2 = \chi_{0.05,1}^2 = 3.841\)

  • Keputusan

    Karena 27.931 > 3.841, maka tolak H0

  • Interpretasi

    Dengan taraf nyata 5%, menolak H0 atau dapat dikatakan bahwa ada interaksi antara jenis kelamin dan pendapat tentang hukuman mati. Dengan kata lain, model yang terbentuk adalah model parameter \(\lambda_{ij}^{XY}\)

  1. Hitung Selisih Deviance
# Deviance of Model
Deviance.model <- model_conditional_Z$deviance - model_homogenous$deviance   # model_conditional_Z: conditional on Z, model_homogenous: homogenous
Deviance.model
## [1] 27.93095
  1. Hitung Derajat Bebas
# Chi Square tabel dengan alpha = 0.05
derajat.bebas <- (3 - 2)
derajat.bebas
## [1] 1
  1. Chi-Square Tabel
chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 3.841459
  1. Keputusan Uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Tolak"

Interpretasi Karena nilai Deviance.model = 27.93 jauh lebih besar dari nilai kritis chi-square tabel = 3.84 (dengan df = 1, alpha = 0.05), maka keputusan uji adalah “Tolak”.

Pada taraf nyata 5%, terdapat cukup bukti untuk menolak H0, atau dengan kata lain ada interaksi antara jenis kelamin (X) dan pendapat tentang hukuman mati (Y). Model terbaik yang terbentuk adalah model dengan menyertakan parameter interaksi (conditional on Z).

Pemilihan Model Terbaik

Ringkasan Model Log Linear

Model Parameter Deviance Jumlah Parameter df AIC
Saturated λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ 0.000 12 0 100.14
Homogenous λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ 1.798 10 2 97.934
Conditional on X λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ 3.9303 8 4 96.067
Conditional on Y λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ 2.9203 8 4 95.057
Conditional on Z λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ 29.729 9 3 123.87

Ringkasan Pengujian Interaksi 3 Arah dan 2 Arah

Interaksi Pengujian \(\Delta\)Deviance \(\Delta\)df Chi-Square Tabel Keputusan Keterangan
XYZ Saturated vs Homogenous 1.798 2 5.991 Terima H0 tidak ada interaksi
YZ Conditional on X vs Homogenous 2.1323 2 5.991 Terima H0 tidak ada interaksi
XZ Conditional on Y vs Homogenous 1.1223 2 5.991 Terima H0 tidak ada interaksi
XY Conditional on Z vs Homogenous 27.931 1 3.841 Tolak H0 ada interaksi

Kesimpulan Pemilihan Model Terbaik

Dari hasil di atas diketahui bahwa asosiasi yang nyata hanya terdapat antara jenis kelamin dan pendapat mengenai hukuman mati. Sehingga, model terbaik adalah:

\[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} \]

Model terbaik adalah model log-linear tanpa interaksi tiga arah dan hanya memuat interaksi dua arah antara jenis kelamin (X) dan sikap terhadap hukuman mati (Y).

Model Terbaik

Model terbaik dipilih berdasarkan pengujian interaksi yang signifikan, yaitu hanya interaksi dua arah antara jenis kelamin (X) dan sikap terhadap hukuman mati (Y):

\[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} \]

# Model Terbaik
bestmodel <- glm(counts ~ x.sex + y.fav + z.fund +
                 x.sex*y.fav,
                 family = poisson(link = "log"))
summary(bestmodel)
## 
## Call:
## glm(formula = counts ~ x.sex + y.fav + z.fund + x.sex * y.fav, 
##     family = poisson(link = "log"))
## 
## Coefficients:
##                   Estimate Std. Error z value Pr(>|z|)    
## (Intercept)        4.26518    0.07794  54.721  < 2e-16 ***
## x.sex1M           -0.59345    0.10645  -5.575 2.48e-08 ***
## y.fav1fav          0.48302    0.08075   5.982 2.20e-09 ***
## z.fund1fund        0.01986    0.07533   0.264    0.792    
## z.fund2mod         0.38130    0.06944   5.491 4.00e-08 ***
## x.sex1M:y.fav1fav  0.65845    0.12708   5.181 2.20e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 245.3612  on 11  degrees of freedom
## Residual deviance:   4.6532  on  6  degrees of freedom
## AIC: 92.79
## 
## Number of Fisher Scoring iterations: 4
Model AIC
Saturated 100.14
Homogenous 97.934
Conditional on X 96.067
Conditional on Y 95.057
Conditional on Z 123.87
Model Terbaik 92.79

Dari hasil di atas, terlihat bahwa model terbaik memiliki AIC yang lebih rendah dibandingkan dengan model lainnya seperti saturated, homogennous, dan conditional model.

Interpretasi Koefisien Model Terbaik

# Interpretasi koefisien model terbaik
data.frame(
  koef = bestmodel$coefficients,
  exp_koef = exp(bestmodel$coefficients)
)
Koefisien Nilai Koefisien Exp(Nilai Koefisien)
Intercept 4.265 71.1776
x.sex1M -0.5934 0.5524
y.fav1fav 0.483 1.621
z.fund1fund 0.0198 1.02
z.fund2mod 0.3813 1.464
x.sex1M:y.fav1fav 0.685 1.9318

Interpretasi Koefisien Model Terbaik

Nilai Dugaan Model Terbaik

Secara manual, nilai fitted value diperoleh dengan cara sebagai berikut:

\[\begin{align*} \\ \hat{\mu}_{111} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{1fav} + \lambda^{z}_{fund} + \lambda^{xy}_{1m,1fav} \right) \\ \\ &= \exp\left( 4.265 - 0.593 + 0.483 + 0.01986 + 0.658 \right) \\ \\ &= \exp(4.833) = 125.595 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{112} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{1fav} + \lambda^{z}_{2mod} + \lambda^{xy}_{1m,1fav} \right) \\ \\ &= \exp\left( 4.265 - 0.593 + 0.483 + 0.381 + 0.658 \right) \\ \\ &= \exp(5.195) = 180.279 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{113} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{1fav} + \lambda^{z}_{3lib} + \lambda^{xy}_{1m,1fav} \right) \\ \\ &= \exp\left( 4.265 - 0.593 + 0.483 + 0 + 0.658 \right) \\ \\ &= \exp(4.813) = 123.126 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{121} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{2opp} + \lambda^{z}_{fund} + \lambda^{xy}_{1m,2opp} \right) \\ \\ &= \exp\left( 4.265 - 0.593 + 0 + 0.01986 + 0 \right) \\ \\ &= \exp(3.692) = 40.109 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{122} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{2opp} + \lambda^{z}_{2mod} + \lambda^{xy}_{1m,2opp} \right) \\ \\ &= \exp\left( 4.265 - 0.593 + 0 + 0.381 + 0 \right) \\ \\ &= \exp(4.053) = 57.572 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{123} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{2opp} + \lambda^{z}_{3lib} + \lambda^{xy}_{1m,2opp} \right) \\ \\ &= \exp\left( 4.265 - 0.593 + 0 + 0 + 0 \right) \\ \\ &= \exp(3.672) = 39.320 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{211} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{1fav} + \lambda^{z}_{fund} + \lambda^{xy}_{2f,1fav} \right) \\ \\ &= \exp\left( 4.265 + 0 + 0.483 + 0.01986 + 0 \right) \\ \\ &= \exp(4.768) = 117.691 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{212} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{1fav} + \lambda^{z}_{2mod} + \lambda^{xy}_{2f,1fav} \right) \\ \\ &= \exp\left( 4.265 + 0 + 0.483 + 0.381 + 0 \right) \\ \\ &= \exp(5.1295) = 168.933 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{213} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{1fav} + \lambda^{z}_{3lib} + \lambda^{xy}_{2f,1fav} \right) \\ \\ &= \exp\left( 4.265 + 0 + 0.483 + 0 + 0 \right) \\ \\ &= \exp(4.748) = 115.377 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{221} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{2opp} + \lambda^{z}_{fund} + \lambda^{xy}_{2f,2opp} \right) \\ \\ &= \exp\left( 4.265 + 0 + 0 + 0.01986 + 0 \right) \\ \\ &= \exp(4.285) = 72.605 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{222} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{2opp} + \lambda^{z}_{2mod} + \lambda^{xy}_{2f,2opp} \right) \\ \\ &= \exp\left( 4.265 + 0 + 0 + 0.381 + 0 \right) \\ \\ &= \exp(4.646) = 104.217 \\ \end{align*}\] \[\begin{align*} \\ \hat{\mu}_{223} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{2opp} + \lambda^{z}_{3lib} + \lambda^{xy}_{2f,2opp} \right) \\ \\ &= \exp\left( 4.265 + 0 + 0 + 0 + 0 \right) \\ \\ &= \exp(4.265) = 71.178 \\ \end{align*}\]

Keterangan: Nilai \(\hat{\mu}_{ijk}\) akan sama apapun referensi dari kategori peubahnya yang digunakan.

Dengan menggunakan R:

# Fitted values dari model terbaik
data.frame(
  Fund = z.fund,
  sex = x.sex,
  favor = y.fav,
  counts = counts,
  fitted = bestmodel$fitted.values
)