KATA PENGANTAR

Pertama-tama saya panjatkan puji syukur ke hadirat Tuhan Yang Maha Esa karena atas berkat dan rahmat-Nya, e-book ini dapat saya susun dan selesaikan.

Dalam era serba digital ini, Artificial Intelligent dan Machune Learning bukanlah hal yang baru lagi. Tidak dapat dipungkiri bahwa saat ini kita hidup di era dimana kecerdasan buatan tidak lagi susah untuk ditemui. Salah satu bentuk pemanfaatan AI juga tercermin dalam proses penyusunan e-book ini. Penggunaan AI bukan tanpa sebab, namun sebagai bentuk bahwa AI dapat digunakan sebagai alat bantu dalam kegiatan belajar mengajar.

E-book ini membahas tentang analisis data di bidang data kategori. Terdapat beberapa materi tentang data kategori seperti: Tabel Kontingensi, Regresi Logistik, Ukuran Asosiasi, Model Regresi Logistik, dan Model Log-Linear.

Melalui e-book ini, saya harap banyak mahasiswa maupun praktisi yang terbantu dalam menganalisis data kategori. Akhir kata, saya menyadari bahwa e-book ini memiliki keterbatasan dalam penulisannya, sehingga kritik dan saran akan sangat membantu saya dalam mengembangkan hal yang lebih baik.

Pendahuluan

Dalam statistik, data kategori adalah jenis data yang dibagi ke dalam kelompok atau kategori dan bukan berupa angka yang bisa diukur. Data seperti ini sering digunakan di berbagai bidang untuk memahami pola, hubungan, dan tren yang tidak bisa dilihat langsung melalui angka. Contoh dari data kategori yaitu jenis kelamin (laki-laki/perempuan), status pernikahan (menikah/belum menikah), tingkat pendidikan (rendah/menengah/tinggi), serta preferensi konsumen (setuju/netral/tidak setuju). Analisis terhadap data kategori sangat penting untuk mengeksplorasi dan memahami informasi yang bersifat nominal atau ordinal.

Karena berbeda dari data numerik, data kategori perlu dianalisis dengan metode khusus, seperti tabel kontingensi, uji chi-square, regresi logistik, dan pendekatan berbasis probabilitas. Seiring berkembangnya teknologi, analisis data kategori juga makin sering digunakan dalam machine learning dan kecerdasan buatan.

Tujuan Analisis Data Kategori

Analisis Data Kategori memiliki banyak tujuan, diantaranya yaitu:

Melihat Frekuensi Data Tersebar dalam Tiap Kategori

Analisis data kategori memungkinkan peneliti untuk melihat bagaimana distribusi data tiap kategori. Misalkan terdapat responden dengan kategori usia muda, dewasa, dan tua. Peneliti dapat melihat sebaran usia dari responden yang ada.

Mengindentifikasi Pola atau Tren

Dalam analisis data kategori, identifikasi pola bertujuan untuk menemukan kecenderungan perbedaan antar kelompok dalam data. Misalnya, melalui visualisasi seperti diagram batang dan diagram lingkaran, kita dapat melihat kecenderungan responden memilih suatu merk dibandingkan dengan merk lain.

Mendukung dalam Pengambilan Keputusan

Dari tujuan yang telah dijelaskan sebelumnya, peneliti juga dimudahkan dalam pengambilan keputusan karena tren atau polanya sudah diketahui.

Definisi dan Ruang Lingkup Analisis Data Kategori

Analisis data kategorik merupakan suatu teknik statistik yang digunakan untuk menganalisis variabel yang bersifat kategorik, dimana data yang dianalisis bersifat nominal atau ordinal.

Definisi Nominal dan Ordinal

Data yang bersifat nominal tidak memiliki urutan, contohnya seperti jenis kelamin dan status menikah. Sedangkan data yang bersifat ordinal memiliki urutan, contohnya seperti tingkat pendidikan dan preferensi konsumen.

Data Biner dan Multikategori

Dalam analisis kategorik, data juga dapat dibedakan berdasarkan jumlah kategorinya. Data biner berarti data hanya memiliki dua kategori, misalnya ya/tidak. Sedangkan data multikategori berarti data tersebut memiliki lebih dari dua kategori.

Perbedaan Data Kategorik dan Data Numerik

Data Kategorik Data Kuantitatif
Bentuk Kategori Numerik
Sifat Kualitatif Kuantitatif
Skala Pengukuran Skala nominal/ordinal Skala Interval/ratio

Manfaat Analisis Data Kategorik dalam Berbagai Bidang

Analisis data kategori memiliki berbagai manfaat dalam berbagai bidang, diantaranya yaitu:

Bidang Tujuan Penerapan
Kesehatan Mengetahui hubungan antara kondisi pasien dan kategori risiko Melihat hubungan antara pola hidup (Sehat/Tidak Sehat) dengan penderita Diabetes
Sosial Menganalisis opini masyarakat Menganalisis kepuasan terhadap layanan publik (Sangat Puas – Tidak Puas)
Ekonomi Mengkategorikan jenis pekerjaan dan status ekonomi Melihat hubungan antara jenis pekerjaan dan tingkat pengeluaran

Metode dalam Analisis Data Kategori

Dalam analisis data kategori, metode yang digunakan bergantung kepada tujuan penelitian. Ada beberapa metode dalam analisis data kategori, diantaranya yaitu:

Uji Chi Square

Uji chi - square bertujuan untuk menguji hubungan antara dua variabel, contohnya yaitu menguji apakah ada hubungan antara jenis kelamin dan preferensi produk.

Uji Fisher’s Exact

Uji Fisher’s Exact bertujuan sama seperti uji chi - square, hanya saya uji ini digunakan ketika frekuensi data berjumlah kecil (<5).

Regresi Logistik

Regresi logistik pada analisis data kategorik bertujuan untuk memprediksi kemungkinan suatu kategori terjadi berdasarkan variabel prediktor. Misalnya memprediksi apakah seseorang siap menikah berdasarkan tingkat pendapatan, kesiapan mental dan fisik, serta kematangan emosional.

Tabulasi Silang (Cross tabulation)

Tabulasi silang atau cross tabulation memiliki tujuan untuk menampilkan distribusi frekuensi dari dua variabel dalam bentuk tabel. Misalnya menampilkan jumlah laki - laki dan perempuan dalam tiap kelompok penyakit.

Analisis Diskriminan (Untuk Data Kategorik Sebagai Target)

Analisis diskriminan bertujuan untuk mengklasifikasikan objek ke dalam kategori berdasarkan beberapa variabel, misalnya menentukan preferensi produk pelanggan jenis kelamin dan pendapatan.

Analisis data kategorik memiliki banyak metode untuk digunakan sesuai tujuan masing-masing. Oleh karena itu, peneliti harus memiliki tujuan dan pemahaman untuk bisa meneliti data kategorik menggunakan metode yang tepat.

Distribusi Probabilitas dalam Data Kategori

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

Distribusi Bernoulli

Distribusi Bernoulli adalah distribusi probabilitas diskrit yang hanya memiliki dua kemungkinan hasil: sukses (1) atau gagal (0). Distribusi ini digunakan untuk percobaan yang hanya memiliki dua hasil, seperti lemparan koin (kepala atau ekor).

Distribusi Bernoulli memiliki satu parameter, yaitu:

  • \(p\): probabilitas sukses (nilai antara 0 dan 1)

Jika \(X\) adalah peubah acak yang mengikuti distribusi Bernoulli, maka:

\[ X \sim \text{Bernoulli}(p) \]

Rumus Distribusi Bernoulli

Probability Mass Function (PMF) dari distribusi Bernoulli didefinisikan sebagai:

\[ P(X = x) = \begin{cases} p & \text{jika } x = 1 \\\\ 1 - p & \text{jika } x = 0 \end{cases} \]

atau bisa ditulis dengan persamaan:

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

Perhitungan dengan R

set.seed(123) # untuk replikasi hasil
n <- 100
hasil <- rbinom(n, size = 1, prob = 0.8)

# Lihat 10 hasil pertama
head(hasil, 10)
##  [1] 1 1 1 0 0 1 1 0 1 1

Distribusi Binomial

Distribusi Binomial adalah distribusi probabilitas diskrit yang menggambarkan jumlah keberhasilan dalam sejumlah percobaan Bernoulli yang independen dan identik. Setiap percobaan memiliki dua hasil: sukses atau gagal, dan probabilitas sukses tetap di setiap percobaan.

Jika \(X\) adalah peubah acak yang mengikuti distribusi binomial, maka:

\[ X \sim \text{Binomial}(n, p) \]

dengan:

  • \(n\): jumlah percobaan
  • \(p\): probabilitas sukses pada setiap percobaan

Rumus Distribusi Binomial

Probability Mass Functiion (PMF) dari distribusi binomial yaitu:

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

Dengan:

  • \(\binom{n}{k}\) adalah koefisien binomial = \(\dfrac{n!}{k!(n-k)!}\)
  • \(k\): jumlah keberhasilan

Perhitungan di R

set.seed(123) # agar hasil bisa direproduksi
simulasi <- rbinom(1000, size = 10, prob = 0.9)
head(simulasi, 10)
##  [1] 10  8  9  8  7 10  9  8  9  9

Distribusi Multinomial

Distribusi Multinomial adalah perluasan dari distribusi Binomial yang mamiliki lebih dari dua hasil. Distribusi ini digunakan untuk menghitung probabilitas hasil dari sejumlah percobaan independen dimana tiap percobaan memiliki lebih dari dua kemungkinan hasil.

Jika terdapat \(k\) kategori dan setiap percobaan memiliki probabilitas:

\[ p_1, p_2, ..., p_k \quad \text{dengan} \quad \sum_{i=1}^k p_i = 1 \]

dan dilakukan \(n\) percobaan, maka jumlah hasil dari setiap kategori mengikuti distribusi multinomial sebagai berikut:

\[ (X_1, X_2, ..., X_k) \sim \text{Multinomial}(n; p_1, p_2, ..., p_k) \]


Rumus Distribusi Multinomial

Probability Mass Functiion dari distribusi Multinomial yaitu:

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

Dengan syarat:

  • \(\sum_{i=1}^k x_i = n\)
  • \(\sum_{i=1}^k p_i = 1\)

Perhitungan dengan R

set.seed(123)
multinomial <- rmultinom(n = 1, size = 15, prob = c(0.2, 0.8, 0.6))
multinomial
##      [,1]
## [1,]    1
## [2,]    7
## [3,]    7

Distribusi Poisson

Distribusi Poisson adalah distribusi probabilitas diskrit yang menggambarkan jumlah kejadian dalam suatu interval waktu tertentu, dengan asumsi bahwa:

  • Kejadian terjadi secara acak,
  • Tidak saling bergantung satu sama lain,
  • Rata-rata kejadian per interval konstan.

Jika \(X\) adalah peubah acak yang mengikuti distribusi Poisson dengan parameter \(\lambda\), maka:

\[ X \sim \text{Poisson}(\lambda) \]

dimana:

  • \(\lambda\): rata-rata jumlah kejadian dalam suatu interval.

Rumus Distribusi Poisson

Probability Mass Functiion (PMF) dari distribusi Poisson adalah:

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

Sifat-sifat penting:

  • Mean: \(\mu = \lambda\)
  • Varians: \(\sigma^2 = \lambda\)

Perhitungan dengan R

set.seed(123)
poisson <- rpois(15, lambda = 5) 
poisson
##  [1] 4 7 4 8 9 2 5 8 5 5 9 5 6 5 2

Desain Sampling dalam Analisis Data Kategori

Desain sampling adalah tahap awal dalam proses penelitian, terkhusus dalam pengumpulan data. Tujuannya adalah agar data yang dikumpulkan dapat digunakan untuk menarik kesimpulan yang akurat tentang populasi.

Dalam analisis data kategori, desain sampling merujuk pada metode yang digunakan untuk memilih sampel dari populasi agar data yang diperoleh dapat mewakili keseluruhan populasi dengan baik. Desain sampling sangatlah penting karena pemilihan sampel akan memengaruhi validitas hasil penelitian.

Desain sampling dalam analisis data kategori dapat diklasifikasikan menjadi dua pendekatan utama, yaitu prospective sampling dan retrosprective sampling.

Prospective Sampling

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

Eksperimen

Dalam studi eksperimen, peneliti menentukan intervensi, yaitu peserta secara acak dibagi ke dalam kelompok perlakuan (menerima intervensi) dan kelompok kontrol (tidak menerima intervensi atau menerima plasebo). Teknik sampling yang digunakan meliputi:

  • Simple Random Sampling (SRS) : Setiap individu memiliki probabilitas yang sama untuk dipilih.
  • Stratified Random Sampling : Populasi dibagi menjadi strata berdasarkan karakteristik tertentu, lalu sampel diambil secara acak dari tiap strata.
  • Clustered Sampling : Populasi dibagi menjadi cluster, kemudian beberapa cluster dipilih secara acak.

Studi Kohort

Studi kohort adalah jenis penelitian dimana peneliti mengamati peserta tanpa intervensi langsung dimana peserta dikelompokkan berdasarkan paparan alami terhadap faktor risiko (misalnya, merokok atau tidak merokok. Jenis sampling yang umum dalam studi yaitu:

  • Systematic Sampling : Subjek dipilih berdasarkan interval tertentu dari populasi
  • Census Sampling : Seluruh anggota dalam populasi tertentu diikutsertakan dalam penelitian

Retrospective Sampling

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

Studi Kasus-Kontrol

Case-Control Study adalah salah satu jenis studi observasional yang digunakan untuk menyelidiki hubungan antara faktor risiko dan suatu outcome (misalnya, penyakit). Berbeda dengan studi kohort, studi ini dimulai dari outcome dan bekerja mundur untuk menentukan paparan sebelumnya. 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.

Studi Kohort Retrospektif

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

  • Quota Sampling : Sampel dipilih untuk mencerminkan proporsi tertentu yang ingin diteliti dalam populasi.
  • Convenience Sampling : Subjek dipilih berdasarkan ketersediaan data yang sudah ada.

Tabel Kontingensi 2 x 2

Tabel kontingensi 2x2 adalah tabel yang digunakan untuk menyajikan data kategorik dari dua variabel yang masing-masing memiliki dua kategori. Tabel ini berguna untuk melihat hubungan atau asosiasi antara dua variabel, misalnya perhitungan risiko relatif (relative risk), perbedaan risiko (risk difference) dan odds ratio.

Struktur umum Tabel Kontingensi 2x2:

Kategori B1 Kategori B2 Total
Kategori A1 a b a + b
Kategori A2 c d c + d
Total a + c b + d n

Distribusi Peluang dalam Tabel Kontingensi 2x2

Peluang Bersama

Peluang bersama adalah probabilitas 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 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.}} \]

Contoh Penerapannya dalam Studi Kasus

Contoh perhitungan: Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.

Penggunaan Sabuk Pengaman Cedera Fatal Cedera Tidak Fatal Total
Tidak memakai sabuk pengaman 1,601 162,527 164,128
Memakai sabuk pengaman 510 412,368 412,878
Total 2,111 574,895 577,006

Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988

Langkah 1 : Hitung Peluang Bersama

  • \[P(\text{No Seatbelt, Fatal}) = \frac{1601}{577006} = 0.0028\]
  • \[P(\text{No Seatbelt, NonFatal}) = \frac{162527}{577006} = 0.2816\]
  • \[P(\text{Seatbelt, Fatal}) = \frac{510}{577006} = 0.0009\]
  • \[P(\text{Seatbelt, NonFatal}) = \frac{412368}{577006} = 0.7147\]

Langkah 2 : Hitung Peluang Marginal

  • \[P(\text{No Seatbelt}) = \frac{164128}{577006} = 0.2844\]
  • \[P(\text{Seatbelt}) = \frac{412878}{577006} = 0.7156\]
  • \[P(\text{Fatal}) = \frac{2111}{577006} = 0.0037\]
  • \[P(\text{NonFatal}) = \frac{574895}{577006} = 0.9963\]

Langkah 3 : Hitung Peluang Bersyarat

  • \[ P(\text{Fatal} \mid \text{No Seatbelt}) = \frac{1601}{164128} = 0.0098 \]

  • \[ P(\text{Fatal} \mid \text{Seatbelt}) = \frac{510}{412878} = 0.0012 \]

  • \[ P(\text{Seatbelt} \mid \text{Fatal}) = \frac{510}{2111} = 0.2416 \]

  • \[ P(\text{Seatbelt} \mid \text{NonFatal}) = \frac{412368}{574895} = 0.7173 \]

Perhitungan dengan R

#Data
data <- matrix(c(1601,162527,510,412368), nrow = 2, byrow = TRUE)
colnames(data) <- c("Fatal", "NonFatal")
rownames(data) <- c("No Seatbelt","Seatbelt")
n <- sum(data)

#Peluang Bersama
p_bersama <- data/n

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

#Peluang Bersyarat
# Peluang bersyarat: P(Fatal | Seatbelt) dan P(Fatal | No Seatbelt)
p_fatal_given_no_seatbelt <- data[1,1] / sum(data[1,])
p_fatal_given_seatbelt <- data[2,1] / sum(data[2,])

# Peluang bersyarat: P(Seatbelt | Fatal) dan P(Seatbelt | NonFatal)
total_fatal <- sum(data[,1])
total_nonfatal <- sum(data[,2])

p_seatbelt_given_fatal <- data[2,1] / total_fatal
p_seatbelt_given_nonfatal <- data[2,2] / total_nonfatal

# Output
p_fatal_given_no_seatbelt
## [1] 0.009754582
p_fatal_given_seatbelt
## [1] 0.001235232
p_seatbelt_given_fatal
## [1] 0.2415917
p_seatbelt_given_nonfatal
## [1] 0.7172927
#Hasil
list(Peluang_Bersama = p_bersama, Peluang_Marginal_Baris = p_marginal_baris, Peluang_Marginal_Kolumn = p_marginal_kolumn)
## $Peluang_Bersama
##                   Fatal  NonFatal
## No Seatbelt 0.002774668 0.2816730
## Seatbelt    0.000883873 0.7146685
## 
## $Peluang_Marginal_Baris
## No Seatbelt    Seatbelt 
##   0.2844476   0.7155524 
## 
## $Peluang_Marginal_Kolumn
##       Fatal    NonFatal 
## 0.003658541 0.996341459

Intepretasi :

Dikarenakan \[ \mathit{P}(\text{Fatal} \mid \text{No Seatbelt}) > \mathit{P}(\text{Fatal} \mid \text{Seatbelt}) \] maka dapat disimpulkan bahwa penggunaan sabuk pengaman dapat menurunkan risiko fatalitas dalam kecelakaan.

Ukuran Asosiasi dalam Data Kategori 2 x 2

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

Risk Difference (RD)

Risk Difference 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 𝐷>0, maka risiko kejadian lebih tinggi di Kelompok 1 dibandingkan Kelompok 2

  • Jika 𝐷<0, maka risiko kejadian lebih rendah di Kelompok 1 dibandingkan Kelompok 2

  • Jika \(\mathit{RD} = 0\), maka tidak ada perbedaan risiko antara dua kelompok.

Relative Risk (RR)

Relative Risk membandingkan kemungkinan terjadinya suatu kejadian pada kelompok yang terpapar terhadap kelompok yang tidak terpapar.

Rumus:

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

  • Jika >1, maka kejadian lebih sering terjadi di Kelompok 1 dibandingkan Kelompok 2.

  • Jika <1, maka kejadian lebih jarang terjadi di Kelompok 1 dibandingkan Kelompok 2.

  • Jika \(\mathit{RR} = 1\), maka tidak ada perbedaan risiko antara dua kelompok.

Odds Ratio (OR)

Odds Ratio digunakan untuk membandingkan peluang terjadinya suatu kejadian antara dua kelompok.

Rumus:

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

Dengan

  • Jika 𝑂 >1, maka peluang kejadian lebih besar di Kelompok 1 dibandingkan Kelompok 2.

  • Jika 𝑂 <1, maka peluang kejadian lebih kecil di Kelompok 1 dibandingkan Kelompok 2.

  • Jika \(\mathit{OR} = 1\), maka tidak ada perbedaan peluang kejadian antara dua kelompok.

Contoh Penerapan dalam Studi Kasus

Contoh perhitungan: Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.

Penggunaan Sabuk Pengaman Cedera Fatal Cedera Tidak Fatal Total
Tidak memakai sabuk pengaman 1,601 162,527 164,128
Memakai sabuk pengaman 510 412,368 412,878
Total 2,111 574,895 577,006

Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988

Contoh Perhitungan Manual

  • Risk Difference \[ RD = \frac{a}{a + b} - \frac{c}{c + d} = \frac{1601}{164128} - \frac{510}{412878} = 0.0098 - 0.0012 = 0.0086 \]

  • Relative Risk \[ \frac{P(\text{Fatal} \mid \text{No Seatbelt})}{P(\text{Fatal} \mid \text{Seatbelt})} = \frac{\frac{1601}{164128}}{\frac{510}{412878}} = \frac{0.009759}{0.001235} \approx 7.897 \]

  • Odds Ratio \[ OR = \frac{a \cdot d}{b \cdot c} = \frac{1601 \cdot 412368}{162527 \cdot 510} = \frac{660201168}{82988770} = 7.96 \]

Contoh Perhitungan dengan R

# Data
a <- 1601
b <- 162527
c <- 510
d <- 412368

# 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.00851935
## 
## $Relative_Risk
## [1] 7.896965
## 
## $Odds_Ratio
## [1] 7.964905

Interpretasi :

  • RD = 0.0086 → Perbedaan risiko kecelakaan fatal sebesar 0.86% lebih tinggi pada pengendara yang tidak pakai sabuk.

  • RR = 7.90 → Pengendara yang tidak memakai sabuk memiliki 8.17 kali risiko lebih besar dalam mengalami kecelakaan fatal.

  • OR = 7.96 → Peluang kecelakaan fatal 7.96 kali lebih besar pada pengendara yang tidak memakai sabuk.

Inferensi Tabel Kontingensi Dua Arah

Inferensi pada statistik merujuk pada proses penarikan kesimpulan tentang suatu populasi berdasarkan data dari sampel. Dalam konteks tabel kontingensi dua arah, inferensi digunakan untuk mengevaluasi dan memahami hubungan antara dua variabel kategorikal yang ditampilkan dalam bentuk tabel matriks.

Tabel kontingensi sendiri menyajikan distribusi frekuensi dari kombinasi kategori dua variabel. Tujuan utama dari inferensi pada tabel kontingensi adalah untuk menganalisis keterkaitan antarvariabel dengan melakukan proses estimasi dan pengujian hipotesis. Secara umum, proses inferensi dalam tabel kontingensi dua arah mencakup dua pendekatan utama yaitu Estimasi dan Pengujian

Estimasi

Estimasi berfungsi untuk memperkirakan parameter populasi berdasarkan sampel. Estimasi dibagi menjadi dua, yaitu:

Estimasi Titik

Estimasi titik digunakan untuk menentukan satu nilai spesifik sebagai perkiraan terbaik dari parameter populasi \[ \hat{p} = \frac{x}{n} \] Dengan:

  • \(\hat{p}\): Estimasi titik proporsi
  • \(x\): Jumlah individu dalam kategori tertentu
  • \(n\): Total jumlah individu dalam sampel

Estimasi Interval

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

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

  • \(Z_{\alpha/2}\): Nilai dari distribusi normal standar untuk tingkat kepercayaan tertentu
  • \(\hat{p}\): Estimasi titik proporsi
  • \(n\): Ukuran sampel

Uji Hipotesis

Uji Proporsi

Uji ini digunakan untuk membandingkan proporsi kejadian antara dua kelompok dalam tabel kontingensi, terutama untuk menentukan apakah terdapat perbedaan yang signifikan dalam proporsi kejadian antara dua kelompok yang berbeda. Tabel kontingensi 2×2 memiliki struktur sebagai berikut:

Kejadian (+) Tidak Kejadian (-) Total (n)
Group 1 n11 n12 n1.
Group 2 n21 n22 n2.
Total n.1 n.2 n

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 diketahui, kita bisa 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.

Perhitungan Manual

Contoh perhitungan: Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.

Penggunaan Sabuk Pengaman Cedera Fatal Cedera Tidak Fatal Total
Tidak memakai sabuk pengaman 1,601 162,527 164,128
Memakai sabuk pengaman 510 412,368 412,878
Total 2,111 574,895 577,006

Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988

  1. Proporsi Masing-Masing Kelompok

Proporsi pada No Seatbelt:

\[ \hat{p}_1 = \frac{1601}{164128} \approx 0.0097 \]

Proporsi pada Seatbelt:

\[ \hat{p}_2 = \frac{510}{412878} \approx 0.0012 \]


  1. Proporsi Gabungan

\[ \hat{p} = \frac{1601+510}{162527+412368} = \frac{2111}{574895} \approx 0.0036 \]


  1. Statistik Uji Z

Rumus:

\[ Z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{\hat{p}(1 - \hat{p})\left(\frac{1}{n_1} + \frac{1}{n_2} \right)}} \]

\[ Z = \frac{0.0097 - 0.0012}{\sqrt{0.0036(1 - 0.0036)\left(\frac{1}{164128} + \frac{1}{412878} \right)}} \approx 48.35712 \]

Perhitungan dengan R

x1 <- 1601
n1 <- 164128
x2 <- 510
n2 <- 412878

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

#Proporsi gabungan
p_pool <- (x1 + x2) / (n1 + n2)
p_pool
## [1] 0.003658541
#Statistik Uji Z
Z <- (p1 - p2) / sqrt(p_pool * (1 - p_pool) * (1/n1 + 1/n2))
Z
## [1] 48.35712
# Tampilkan hasil
p1
## [1] 0.009754582
p2
## [1] 0.001235232
p_pool
## [1] 0.003658541
Z
## [1] 48.35712

Interpretasi : Jika kita menggunakan taraf signifikansi sebesar 5%, kita menerima H0 karena nilai |Z| yaitu 48.35712 > dari 1.96 yang berarti bahwa ada 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.

Tiga ukuran utama dalam uji asosiasi adalah:

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

Perbandingan antara risiko kejadian di kelompok 1 dibanding kelompok 2.

Rumus:

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

Standar Error untuk 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

Perbandingan peluang kejadian antara dua kelompok

Rumus: \[ OR = \frac{a \cdot d}{b \cdot c} \] \[ 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))} \]

Contoh Perhitungan

Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.

Penggunaan Sabuk Pengaman Cedera Fatal Cedera Tidak Fatal Total
Tidak memakai sabuk pengaman 1,601 162,527 164,128
Memakai sabuk pengaman 510 412,368 412,878
Total 2,111 574,895 577,006

Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988

  1. Risk Difference (RD) \[ \hat{p}_1 = \frac{1601}{164128} \approx 0.00976, \quad \hat{p}_2 = \frac{510}{412878} \approx 0.00124 \]

\[ RD = \hat{p}_1 - \hat{p}_2 = 0.00976 - 0.00124 = 0.00852 \]

\[ SE(RD) = \sqrt{ \frac{0.00976(1 - 0.00976)}{164128} + \frac{0.00124(1 - 0.00124)}{412878} } \] \[ = \sqrt{5.89 \times 10^{-8} + 3.01 \times 10^{-9}} = \sqrt{6.19 \times 10^{-8}} = 0.000249 \]

\[ Z = \frac{0.00852}{0.000249} \approx 34.22 \]

  1. Relative Risk (RR)

\[ RR = \frac{0.00976}{0.00124} = 7.87 \]

Standard Error log(RR)

\[ SE(\log(RR)) = \sqrt{ \left( \frac{1}{1601} - \frac{1}{164128} \right) + \left( \frac{1}{510} - \frac{1}{412878} \right) } \]

\[ = \sqrt{(0.0006246 - 0.0000061) + (0.0019608 - 0.0000024)} = \sqrt{0.0006185 + 0.0019584} \] \[ = \sqrt{0.0025769} = 0.05076 \] Statistik uji Z untuk RR

\[ Z_{RR} = \frac{\log(7.87)}{0.05076} = \frac{2.06}{0.05076} \approx 40.58 \]

  1. Odds Ratio

\[ OR = \frac{a \cdot d}{b \cdot c} = \frac{1601 \cdot 412368}{162527 \cdot 510} = \frac{660242368}{82888770} = 7.97 \] Standard Error log(OR)

\[ SE(\log(OR)) = \sqrt{ \frac{1}{1601} + \frac{1}{162527} + \frac{1}{510} + \frac{1}{412368} } = \sqrt{0.0006246 + 0.0000062 + 0.0019608 + 0.0000024} = \sqrt{0.0025939} = 0.05093 \] Statistik uji Z untuk OR

\[ Z_{OR} = \frac{\log(7.97)}{0.05093} = \frac{2.076}{0.05093} \approx 40.77 \]

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

di mana:

  • \(O\) adalah nilai observasi dalam tabel kontingensi.
  • \(E\) adalah nilai yang diharapkan, dihitung sebagai:

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

  • \(R_i\) = total baris ke-\(i\),
  • \(C_j\) = total kolom ke-\(j\),
  • \(N\) = total sampel.
Contoh perhitungan Chi Square

Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.

Penggunaan Sabuk Pengaman Cedera Fatal Cedera Tidak Fatal Total
Tidak memakai sabuk pengaman 1,601 162,527 164,128
Memakai sabuk pengaman 510 412,368 412,878
Total 2,111 574,895 577,006

Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988

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

dimana:

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

Menghitung nilai yang diharapkan:

Gunakan rumus:

\[ E_{ij} = \frac{\text{Total baris } i \times \text{Total kolom } j}{\text{Total keseluruhan}} \]

Untuk Tidak Pakai Sabuk, Cedera Fatal (\(E_{11}\)):

\[ E_{11} = \frac{164128 \times 2111}{577006} = \frac{346555408}{577006} \approx 600.2 \]

Untuk Tidak Pakai Sabuk, Tidak Fatal (\(E_{12}\)):

\[ E_{12} = \frac{164128 \times 574895}{577006} \approx 163527.8 \]

Untuk Pakai Sabuk, Cedera Fatal (\(E_{21}\)):

\[ E_{21} = \frac{412878 \times 2111}{577006} \approx 1510.8 \]

Untuk Pakai Sabuk, Tidak Fatal (\(E_{22}\)):

\[ E_{22} = \frac{412878 \times 574895}{577006} \approx 411367.2 \]

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

Gunakan rumus:

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

Untuk Tidak Pakai Sabuk, Fatal:

\[ \frac{(1601 - 600.2)^2}{600.2} = \frac{1000800.64}{600.2} \approx 1667.9 \]

Untuk Tidak Pakai Sabuk, Tidak Fatal:

\[ \frac{(162527 - 163527.8)^2}{163527.8} = \frac{1008010.24}{163527.8} \approx 6.2 \]

Untuk Pakai Sabuk, Fatal:

\[ \frac{(510 - 1510.8)^2}{1510.8} = \frac{1008010.24}{1510.8} \approx 667.2 \]

Untuk Pakai Sabuk, Tidak Fatal:

\[ \frac{(412368 - 411367.2)^2}{411367.2} = \frac{1000800.64}{411367.2} \approx 2.4 \]

Total Chi-Square:

\[ \chi^2 = 1667.9 + 6.2 + 667.2 + 2.4 = 2343.7 \]

Perhitungan menggunakan R

# Membuat data tabel 2x2
data <- matrix(c(1601, 162527, 510, 412368), nrow = 2, byrow = TRUE)

# Menambahkan nama baris dan kolom
rownames(data) <- c("Tidak Memakai Sabuk", "Memakai Sabuk")
colnames(data) <- c("Cedera Fatal", "Cedera Tidak Fatal")

# Menampilkan data tabel
data
##                     Cedera Fatal Cedera Tidak Fatal
## Tidak Memakai Sabuk         1601             162527
## Memakai Sabuk                510             412368
# Menghitung uji 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 = 2336.1, df = 1, p-value < 2.2e-16

Interpretasi: Dengan taraf signifikansi 5%, didapatkan p-value < 0.05, maka kita menolak H0 (hipotesis nol). Artinya, ada hubungan yang signifikan secara statistik antara penggunaan sabuk pengaman dan tingkat fatalitas kecelakaan

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×Jsignifikan, maka partisi Chi-Square memungkinkan kita untuk menguraikan efek hubungan dalam subkelompok yang lebih kecil.

Contoh Perhitungan

Contoh Kasus:

Identifikasi Partai Berdasarkan Ras (GSS 1991)

Ras Demokrat Independen Republik
Black 103 15 11
White 341 105 405
Total 444 120 416

Sumber: 1991 General Social Survey, National Opinion Research Center

Perhitungan Manual Langkah 1: Menghitung Frekuensi yang Diharapkan

Frekuensi yang diharapkan dihitung menggunakan rumus:

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

\[ E_{1,1} = \frac{129 \times 444}{980} = 58.48 \]

\[ E_{1,2} = \frac{129 \times 120}{980} = 15.79 \]

\[ E_{1,3} = \frac{129 \times 416}{980} = 54.73 \]

\[ E_{2,1} = \frac{851 \times 444}{980} = 385.52 \]

\[ E_{2,2} = \frac{851 \times 120}{980} = 104.20 \]

\[ E_{2,3} = \frac{851 \times 416}{980} = 361.27 \]

Langkah 2: Menghitung Statistik Chi-Square

Sekarang kita dapat menghitung statistik chi-square \[\chi^2\] dengan rumus:

\[ \chi^2 = \sum \frac{(O_{ij} - E_{ij})^2}{E_{ij}} \]

Dimana \(O_{ij}\) adalah frekuensi yang diamati dan \(E_{ij}\) adalah frekuensi yang diharapkan.

Dimana \(O_{ij}\) adalah frekuensi yang diamati dan \(E_{ij}\) adalah frekuensi yang diharapkan. Berdasarkan data yang ada, kita menghitung chi-square untuk setiap sel:

\[ \chi^2 = \frac{(103 - 58.48)^2}{58.48} + \frac{(15 - 15.79)^2}{15.79} + \frac{(11 - 54.73)^2}{54.73}+ \] \[ \frac{(341 - 385.52)^2}{385.52} + \frac{(105 - 104.20)^2}{104.20} + \frac{(405 - 361.27)^2}{361.27} \]

\[ = 33.89 + 0.0395 + 34.95 + 5.14 + 0.0061 + 5.29 = 79.32 \]

\[ \chi^2 \approx 79.32 \]

\[ \chi^2 = 33.89 + 0.0395 + 34.95 + 5.14 + 0.0061 + 5.29 = 79.32 \]

Jadi, berdasarkan perhitungan manual:

\[ \chi^2 \approx 79.32 \]

Perhitungan Dengan R:

# Membuat matriks data
data_matrix <- matrix(c(103, 15, 11,
                        341, 105, 405), 
                      nrow = 2, byrow = TRUE)

colnames(data_matrix) <- c("Demokrat", "Independen", "Republik")
rownames(data_matrix) <- c("Black", "White")

# Melakukan uji chi-square
chi_square_test <- chisq.test(data_matrix)

# Menampilkan hasil uji
chi_square_test
## 
##  Pearson's Chi-squared test
## 
## data:  data_matrix
## X-squared = 79.431, df = 2, p-value < 2.2e-16

Interpretasi: Karena nilai \(\chi^2\) sangat besar dan jauh melebihi nilai kritis untuk df = 2 pada taraf signifikansi 0.05, maka kita menolak \(H_0\). Ini menunjukkan bahwa ada hubungan signifikan antara ras dan identifikasi partai.

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

Dimana:

- \(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

Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.

Penggunaan Sabuk Pengaman Cedera Fatal Cedera Tidak Fatal Total
Tidak memakai sabuk pengaman 1,601 162,527 164,128
Memakai sabuk pengaman 510 412,368 412,878
Total 2,111 574,895 577,006

Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988

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{164128 \times 2111}{577006} = 599.78 \\ \hat{\mu}_{12} &= \frac{164128 \times 574895}{577006} = 163528.22 \\ \hat{\mu}_{21} &= \frac{412878 \times 2111}{577006} = 1511.22 \\ \hat{\mu}_{22} &= \frac{412878 \times 574895}{577006} = 411366.78 \\ \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 \bigg[ 1601 \ln \left( \frac{1601}{599.78} \right) + 162527 \ln \left( \frac{162527}{163528.22} \right) + \\ &\quad 510 \ln \left( \frac{510}{1511.22} \right) + 412368 \ln \left( \frac{412368}{411366.78} \right) \bigg] \end{aligned} \]

Setelah dihitung:

\[ \begin{aligned} G^2 &= 2 \times (1693.26 + (-6.14) + (-1037.45) + 2.43) \\ &= 2 \times 652.10 = 1304.20 \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.

Perhitungan Dengan R

# Data Observasi
observed <- matrix(c(1601, 162527,
                     510, 412368), 
                   nrow = 2, byrow = TRUE)

rownames(observed) <- c("Tidak Pakai Sabuk", "Pakai Sabuk")
colnames(observed) <- c("Cedera Fatal", "Tidak Fatal")

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

# Frekuensi harapan
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_critical <- qchisq(0.95, df)

# Output hasil
cat("Pearson's Chi Square =", round(G2, 2), "\n")
## Pearson's Chi Square = 2041.16
cat("Derajat bebas =", df, "\n")
## Derajat bebas = 1
cat("Nilai kritis chi-square (0.05) =", round(chi_critical, 3), "\n")
## Nilai kritis chi-square (0.05) = 3.841
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.

Interpretasi: Karena \(G^2 = 1304.20 > 3.841\), maka kita menolak hipotesis nol, artinya terdapat bukti yang cukup bahwa penggunaan sabuk pengaman berpengaruh signifikan terhadap tingkat fatalitas kecelakaan.

Tabel Kontingensi Tiga Arah

Tabel kontingensi tiga arah adalah perpanjangan dari tabel kontingensi dua arah yang digunakan untuk menganalisis hubungan antara tiga variabel kategori secara simultan. Dalam banyak situasi, hubungan antara dua variabel (misalnya X dan Y) dapat dipengaruhi oleh variabel ketiga Z, yang disebut sebagai variabel kontrol. Tabel ini sering digunakan dalam analisis data ketika ada potensi faktor pengganggu yang dapat mempengaruhi hubungan antara dua variabel utama. Misalnya:

  • Dalam studi epidemiologi, untuk meneliti hubungan antara merokok (X) dan kanker paru-paru (Y) dengan mengendalikan efek usia (Z).
  • Dalam penelitian sosial, untuk mengevaluasi hubungan antara ras tersangka (X) dan keputusan hukuman mati (Y) dengan mempertimbangkan kas korban (Z).

Tabel kontingensi tiga arah dapat dibagi menjadi dua jenis:

  1. Tabel Parsial: Tabel yang menyajikan hubungan antara X dan Y pada setiap kategori Z. Ini memungkinkan analisis hubungan bersyarat, dimana efek Z dikendalikan.

  2. Tabel Marginal: Tabel yang diperoleh dengan mengabaikan Z, yaitu dengan menjumlahkan semua kategori Z. Ini memberikan gambaran umum hubungan antara X dan Y tanpa mempertimbangkan Z, yang bisa menyebabkan distorsi interpretasi, seperti dalam Simpson’s Paradox

Kegunaan tabel marginal adalah untuk melihat pola asosiasi secara agregat, tetapi sering kali mengabaikan efek kovariat yang dapat memberikan pemahaman lebih mendalam. Oleh karena itu, analisis yang mempertimbangkan tabel parsial biasanya lebih ditentukan untuk mendapatkan kesimpulan yang lebih akurat dalam penelitian ilmiah.

Tabel Parsial dan Marginal

Tabel parsial adalah tabel yang mengelompokkan𝑋dan𝑌berdasarkan setiap level 𝑍, sedangkan tabel marginal adalah tabel yang mengabaikan𝑍,dengan menjumlahkan data dari semua level𝑍.

Tabel berikut menunjukkan data mengenai efektivitas obat AZT (zidovudine) dalam memperlambat perkembangan gejala AIDS. Studi ini melibatkan 338 veteran yang rasnya dibagi menjadi dua yaitu ras kulit putih dan ras kulit hitam yang telah terinfeksi virus AIDS dan mulai menunjukkan tanda-tanda penurunan sistem kekebalan tubuh.

Race AZT Use Yes (Symptoms) No (Symptoms)
White Yes 14 93
White No 32 81
Black Yes 11 52
Black No 12 43

Sumber: New York Times, Feb 15, 1991.

Tabel Frekuensi Parsial untuk Z (Race) = White

AZT Use Yes (Symptoms) No (Symptoms)
Yes 14 93
No 32 81

Tabel Frekuensi Parsial untuk Z (Race) = White menunjukkan hubungan antara penggunaan AZT dan perkembangan gejala AIDS pada kelompok ras kulit putih.

Tabel Frekuensi Parsial untuk Z (Race) = Black

AZT Use Yes (Symptoms) No (Symptoms)
Yes 11 52
No 12 43

Tabel Frekuensi Parsial untuk Z (Race) = Black menunjukkan hubungan antara penggunaan AZT dan perkembangan gejala AIDS pada kelompok ras kulit hitam.

Dengan R:

# Membuat data dengan array
data_azt <- array(
  c(14, 93, 32, 81, 11, 52, 12, 43), 
  dim = c(2, 2, 2), 
  dimnames = list(
    AZT_Use = c("Yes", "No"),
    Symptoms = c("Yes", "No"),
    Race = c("White", "Black")
  )
)

# Tampilkan array untuk memastikan datanya
data_azt
## , , Race = White
## 
##        Symptoms
## AZT_Use Yes No
##     Yes  14 32
##     No   93 81
## 
## , , Race = Black
## 
##        Symptoms
## AZT_Use Yes No
##     Yes  11 12
##     No   52 43
# Ekstrak tabel frekuensi parsial berdasarkan ras
freq_parsial_white <- data_azt[, , "White"]
freq_parsial_black <- data_azt[, , "Black"]

# Tampilkan hasil frekuensi parsial
freq_parsial_white
##        Symptoms
## AZT_Use Yes No
##     Yes  14 32
##     No   93 81
freq_parsial_black
##        Symptoms
## AZT_Use Yes No
##     Yes  11 12
##     No   52 43

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 Race dan Symptoms

Race Yes No Total
White 46 174 220
Black 23 95 118
Total 69 269 338

Tabel frekuensi marginal untuk AZT Use dan Symptoms

AZT Use Yes No Total
Yes 25 145 170
No 44 124 168
Total 69 269 338

Dengan R:

# Buat data array 3D: AZT Use x Symptoms x Race
data_azt <- array(
  c(
    14, 32,  # White - AZT Yes
    93, 81,  # White - AZT No
    11, 12,  # Black - AZT Yes
    52, 43   # Black - AZT No
  ),
  dim = c(2, 2, 2),
  dimnames = list(
    AZT_Use = c("Yes", "No"),
    Symptoms = c("Yes", "No"),
    Race = c("White", "Black")
  )
)

# Tampilkan array
data_azt
## , , Race = White
## 
##        Symptoms
## AZT_Use Yes No
##     Yes  14 93
##     No   32 81
## 
## , , Race = Black
## 
##        Symptoms
## AZT_Use Yes No
##     Yes  11 52
##     No   12 43
# Frekuensi marginal untuk AZT Use (menjumlah semua Symptoms & Race)
freq_marginal_azt <- apply(data_azt, 1, sum)

# Frekuensi marginal untuk Symptoms (menjumlah semua AZT Use & Race)
freq_marginal_symptoms <- apply(data_azt, 2, sum)

# Tampilkan hasil
freq_marginal_azt
## Yes  No 
## 170 168
freq_marginal_symptoms
## Yes  No 
##  69 269

Distribusi Peluang

Tabel berikut menunjukkan data mengenai efektivitas obat AZT (zidovudine) dalam memperlambat perkembangan gejala AIDS. Studi ini melibatkan 338 veteran yang rasnya dibagi menjadi dua yaitu ras kulit putih dan ras kulit hitam yang telah terinfeksi virus AIDS dan mulai menunjukkan tanda-tanda penurunan sistem kekebalan tubuh.

Race AZT Use Yes (Symptoms) No (Symptoms)
White Yes 14 93
White No 32 81
Black Yes 11 52
Black No 12 43

Peluang bersama didefinisikan sebagai:

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

Contoh Perhitungan dengan Data diatas

# Load library yang diperlukan
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
# Data dalam format array 3D
data_azt <- array(
  c(
    14, 32,  # White - AZT Yes
    93, 81,  # White - AZT No
    11, 12,  # Black - AZT Yes
    52, 43   # Black - AZT No
  ),
  dim = c(2, 2, 2),
  dimnames = list(
    AZT_Use = c("Yes", "No"),
    Symptoms = c("Yes", "No"),
    Race = c("White", "Black")
  )
)


# Menampilkan tabel array 3D
print(data_azt)
## , , Race = White
## 
##        Symptoms
## AZT_Use Yes No
##     Yes  14 93
##     No   32 81
## 
## , , Race = Black
## 
##        Symptoms
## AZT_Use Yes No
##     Yes  11 52
##     No   12 43
total <- sum(data_azt)
joint_prob <- data_azt/ total
ftable(joint_prob)
##                  Race      White      Black
## AZT_Use Symptoms                           
## Yes     Yes           0.04142012 0.03254438
##         No            0.27514793 0.15384615
## No      Yes           0.09467456 0.03550296
##         No            0.23964497 0.12721893

Interpretasi: Secara umum, data menunjukkan bahwa AZT mungkin efektif dalam menurunkan gejala AIDS, tapi belum bisa disimpulkan secara pasti tanpa analisis statistik lanjutan.

Peluang Marginal

Peluang marginal didefinisikan sebagai berikut

untuk Ras: \[ P(Z) = \frac{n(Z)}{N} \]

marginal_Z <- apply(joint_prob, 3, sum)
marginal_Z
##     White     Black 
## 0.6508876 0.3491124

Interpretasi : Proporsi Ras Putih lebih banyak dibanding ras Hitam.

Untuk AZT use: \[ P(X) = \frac{n(X)}{N} \]

marginal_X <- apply(joint_prob, 1, sum)
marginal_X
##       Yes        No 
## 0.5029586 0.4970414

Intepretasi : Proporsi penggunaan AZT lebih banyak dibandingkan yang tidak memakai AZT

Peluang Bersyarat

Contoh: Mencari peluang seseorang ber-ras tertentu jika ia mendapat AZT dan mengalami gejala AIDS:

\[ P(Race| AZT = Yes, Symptoms = Yes) = \frac{P(Race, AZT = Yes, Symptoms = Yes)}{P(AZT = Yes, Symptoms = Yes)} \]

# Hitung peluang bersyarat P(Race | AZT Use, Symptoms)
p_Race_given_AZT_Symptoms <- prop.table(data_azt, margin = c(1, 2))
p_Race_given_AZT_Symptoms
## , , Race = White
## 
##        Symptoms
## AZT_Use       Yes        No
##     Yes 0.5600000 0.6413793
##     No  0.7272727 0.6532258
## 
## , , Race = Black
## 
##        Symptoms
## AZT_Use       Yes        No
##     Yes 0.4400000 0.3586207
##     No  0.2727273 0.3467742

Interpretasi:

Peluang seseorang ber-ras putih jika ia mendapat AZT dan mengalami gejala AIDS adalah 0.56.

Peluang seseorang ber-ras hitam jika ia mendapat AZT dan mengalami gejala AIDS adalah 0.44.

Ukuran Asosiasi

Perbedaan Peluang (Risk Difference, RD)

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

#Hitung Risk Difference (RD) untuk masing-masing ras

# White
rd_white <- (data_azt["Yes", "Yes", "White"] / sum(data_azt["Yes", , "White"])) -
            (data_azt["No", "Yes", "White"] / sum(data_azt["No", , "White"]))

# Black
rd_black <- (data_azt["Yes", "Yes", "Black"] / sum(data_azt["Yes", , "Black"])) -
            (data_azt["No", "Yes", "Black"] / sum(data_azt["No", , "Black"]))

# Tampilkan hasil
rd_white
## [1] -0.1523447
rd_black
## [1] -0.04357864

Interpretasi: AZT mengurangi risiko berkembangnya gejala AIDS sebesar 15.2% pada ras White. AZT mengurangi risiko sebesar 4.36% pada ras Black, tapi efeknya lebih kecil dibanding pada ras White.

Risiko Relatif (Relative Risk, RR)

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

#Hitung Risk Relative (RR) untuk masing-masing ras

# White
rr_white <- (data_azt["Yes", "Yes", "White"] / sum(data_azt["Yes", , "White"])) /
            (data_azt["No", "Yes", "White"] / sum(data_azt["No", , "White"]))

# Black
rr_black <- (data_azt["Yes", "Yes", "Black"] / sum(data_azt["Yes", , "Black"])) /
            (data_azt["No", "Yes", "Black"] / sum(data_azt["No", , "Black"]))

# Tampilkan hasil
rr_white
## [1] 0.4620327
rr_black
## [1] 0.8002646

Interpretasi: Pada kelompok ras White, risiko berkembangnya gejala AIDS bagi yang menggunakan AZT adalah sekitar 46% dari risiko pada yang tidak menggunakan AZT. Pada kelompok ras Black, risiko berkembangnya gejala AIDS bagi yang menggunakan AZT adalah sekitar 80% dari risiko pada yang tidak menggunakan AZT.

Odds Ratio (OR)

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

# Menghitung Odds Ratio untuk masing masing ras
or_white <- (data_azt["Yes", "Yes", "White"] / data_azt["Yes", "No", "White"]) /
            (data_azt["No", "Yes", "White"] / data_azt["No", "No", "White"])

or_black <- (data_azt["Yes", "Yes", "Black"] / data_azt["Yes", "No", "Black"]) /
            (data_azt["No", "Yes", "Black"] / data_azt["No", "No", "Black"])

or_white
## [1] 0.3810484
or_black
## [1] 0.7580128

Interpretasi: Di antara individu dengan ras putih, peluang (odds) mengalami gejala AIDS pada yang menggunakan AZT adalah sekitar 38% dari odds mengalami gejala pada yang tidak menggunakan AZT Di antara individu dengan ras hitam, peluang (odds) mengalami gejala AIDS pada yang menggunakan AZT adalah sekitar 76% dari odds mengalami gejala pada yang tidak menggunakan AZT

Inferensi Tabel Kontingensi Tiga Arah

Tabel kontingensi tiga arah digunakan untuk menganalisis hubungan antara dua variabel kategorik dengan mempertimbangkan variabel kontrol. Contohnya adalah hubungan antara Penggunaan AZT (X) dan perkembangan gejala AIDS (Y), dengan variabel kontrol ras (Z)

Tabel ini terdiri dari beberapa tabel parsial (2×2) untuk setiap tingkat Z, serta tabel marginal yang mengabaikan Z. Ukuran asosiasi yang digunakan adalah odds ratio.

Jika odds ratio parsial relatif konstan, kita dapat menghitung odds ratio bersama menggunakan estimasi Mantel-Haenszel.

Independensi Bersyarat dalam Tabel Kontingensi Tiga Arah

Independensi bersyarat adalah konsep penting dalam analisis tabel kontingensi tiga arah. Ini merujuk pada kondisi di mana dua variabel, 𝑋 dan 𝑌, independen dalam setiap level variabel ketiga, 𝑍. Pengujian independensi bersyarat dilakukan dengan metode statistik seperti uji Cochran-Mantel-Haenszel (CMH).

Definisi Independensi Bersyarat

Independensi Bersyarat: Dua variabel, 𝑋 dan 𝑌, dikatakan independen bersyarat terhadap variabel ketiga, 𝑍, jika rasio odds mereka dalam setiap strata 𝑍 sama dengan 1.

Pengujian Statistik untuk Independensi Bersyarat

Metode Cochran-Mantel-Haenszel Tujuan Uji CMH Uji Cochran-Mantel-Haenszel (CMH) digunakan untuk menguji hubungan antara dua variabel kategori dengan mempertimbangkan efek dari variabel perancu (confounder).

Statistik uji Cochran-Mantel-Haenszel (CMH) dirumuskan sebagai:

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

Keterangan

  • \(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)} \] Contoh Perhitungan:

Race: White

Symptoms (Yes) Symptoms (No) Total
AZT Yes 14 93 107
AZT No 32 81 113

Race: Black

Symptoms (Yes) Symptoms (No) Total
AZT Yes 11 52 63
AZT No 12 43 55

Perhitungan Manual:

  1. Nilai Harapan \(\mu_{1ik}\)

\[ \mu_{1ik} = \frac{n_{i+k} \cdot n_{+ik}}{n_{++k}} \]

Untuk Ras Putih (k=1):

\[ \mu_{111} = \frac{107 \times 46}{220} = 22.38 \]

Untuk Ras Hitam (k=2):

\[ \mu_{112} = \frac{63 \times 23}{118} = 12.28 \]

  1. Varians

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

Untuk Ras Putih (k=1):

\[ \text{Var}(n_{111}) = \frac{107 \cdot 113 \cdot 46 \cdot 174}{220^2 \cdot 219} \approx 9.41 \]

Untuk Ras Hitam (k=2):

\[ \text{Var}(n_{112}) = \frac{63 \cdot 55 \cdot 23 \cdot 95}{118^2 \cdot 117} \approx 5.25 \]

  1. Statistik CMH

\[ X^2_{CMH} = \frac{(n_{111} - \mu_{111})^2 + (n_{112} - \mu_{112})^2}{\text{Var}(n_{111}) + \text{Var}(n_{112})} \]

\[ X^2_{CMH} = \frac{(14 - 22.38)^2 + (11 - 12.28)^2}{9.41 + 5.25} = \frac{74.16}{14.66} \approx 5.06 \]

Interpretasi: Dengan taraf signifikansi 5%, nilai kritis \(\chi^2\) adalah 3.841. Karena \(X^2_{CMH} = 5.06 > 3.841\), maka hipotesis nol ditolak. Artinya, ada hubungan signifikan antara penggunaan AZT dan perkembangan gejala AIDS setelah mengontrol variabel ras.

Perhitungan dengan R:

# Data dalam format array 3D
data_azt <- array(
  c(
    14, 32,  # White - AZT Yes
    93, 81,  # White - AZT No
    11, 12,  # Black - AZT Yes
    52, 43   # Black - AZT No
  ),
  dim = c(2, 2, 2),
  dimnames = list(
    AZT_Use = c("Yes", "No"),
    Symptoms = c("Yes", "No"),
    Race = c("White", "Black")
  )
)

# Menampilkan tabel
print(data_azt)
## , , Race = White
## 
##        Symptoms
## AZT_Use Yes No
##     Yes  14 93
##     No   32 81
## 
## , , Race = Black
## 
##        Symptoms
## AZT_Use Yes No
##     Yes  11 52
##     No   12 43
cmh_base <- mantelhaen.test(data_azt, correct = FALSE)
cmh_base
## 
##  Mantel-Haenszel chi-squared test without continuity correction
## 
## data:  data_azt
## Mantel-Haenszel X-squared = 6.7624, df = 1, p-value = 0.00931
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
##  0.2818621 0.8414208
## sample estimates:
## common odds ratio 
##         0.4869955

Interpretasi : Dikarenakan p-value < 0.05, maka keputusan uji adalah menolak H₀, yang berarti terdapat hubungan yang signifikan antara penggunaan AZT dengan perkembangan gejala AIDS setelah mengontrol variabel ras. Dengan kata lain, AZT secara signifikan berhubungan dengan timbul atau tidaknya gejala AIDS, terlepas dari perbedaan ras peserta dalam studi ini

Generalized Linear Model (GLM)

Generalized Linear Model (GLM) merupakan perluasan dari model regresi linear klasik. GLM memungkinkan digunakan untuk memodelkan data dimana variabel respons tidak berdistribusi normal dan atau hubungan antara variabel prediktor dengan rata-rata dari variabel respons tidak linear. GLM terdiri dari tiga komponen utama:

  1. Distribusi dari exponential family untuk variabel respons.
  2. Fungsi link yang menghubungkan ekspektasi dari variabel respons ke kombinasi linear dari variabel prediktor.
  3. Fungsi linear prediktor : 𝜂 =𝑋𝛽

Exponential Family

Distribusi termasuk dalam exponential family jika dapat ditulis dalam bentuk:

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

Contoh Distribusi yang Termasuk Exponential Family:

  • Distribusi Normal
  • Distribusi Binomial
  • Distribusi Poisson
  • Distribusi Gamma

Contoh Pembuktian: Distribusi Binomial

Fungsi probabilitas distribusi binomial adalah:

\[ P(Y = y) = \binom{n}{y} \pi^y (1 - \pi)^{n - y} \]

Tuliskan ulang dalam bentuk exponential family:

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

Dengan substitusi:

  • \(\theta = \log \left( \frac{\pi}{1 - \pi} \right)\)
  • \(b(\theta) = -n \log(1 - \pi)\)
  • \(\phi = 1\)

Kesimpulan : Dengan bentuk tersebut, maka distribusi binomial termasuk dalam keluarga exponential family.

Model Regresi Logistik

Regresi Logistik memiliki kemiripan dengan regresi linear, yaitu menggabungkan nilai-nilai input secara linear dengan bobot (koefisien) untuk menghasilkan sebuah prediksi. Bedanya, regresi logistik membatasi hasil prediksinya agar berada dalam bentuk nilai biner 0 atau 1 dengan bantuan fungsi sigmoid sebagai fungsi aktivasi. Fungsi ini mengubah output menjadi nilai antara 0 sampai 1 dan membentuk kurva seperti huruf S.

Model ini digunakan untuk menganalisis hubungan antara satu atau lebih variabel independen dan mengklasifikasikan data ke dalam kelas-kelas yang bersifat diskrit. Artinya, regresi logistik cocok digunakan untuk masalah klasifikasi, terutama yang hanya memiliki dua kemungkinan kategori (dikenal dengan klasifikasi biner).

Contohnya, angka 0 bisa merepresentasikan gagal, sedangkan angka 1 merepresentasikan sukses.

Beberapa contoh penerapan regresi logistik dalam klasifikasi biner antara lain:

-Memprediksi kelulusan siswa: Model regresi logistik dapat digunakan untuk memperkirakan apakah seorang siswa akan lulus atau tidak berdasarkan variabel-variabel seperti kehadiran, nilai tugas, partisipasi kelas, dan nilai ujian.

-Deteksi penipuan transaksi keuangan: Dalam dunia perbankan dan keuangan, regresi logistik dapat membantu mendeteksi apakah sebuah transaksi tergolong normal atau mencurigakan (fraud) dengan menganalisis pola transaksi seperti frekuensi, jumlah uang, lokasi transaksi, dan waktu.

-Memprediksi hasil diagnosis penyakit: Dalam dunia medis, model ini bisa digunakan untuk menentukan apakah pasien menderita suatu penyakit (misalnya diabetes atau kanker) berdasarkan variabel seperti usia, tekanan darah, kadar gula, dan riwayat keluarga.

Keunggulan Utama Regresi Logistik yaitu mudah diterapkan dalam Machine Learning. Regresi logistik sangat cocok digunakan dalam Machine Learning karena proses pelatihan (training) dan pengujian (testing) nya cukup sederhana. Model ini belajar dari pola-pola dalam data input dan menghubungkannya dengan output (label). Karena tidak membutuhkan komputasi tinggi, regresi logistik tergolong mudah untuk diimplementasikan, dipahami, dan dilatih dibanding metode lain dalam Machine Learning.

Cocok untuk data yang bisa dipisahkan secara linear Jika dua kelas data bisa dipisahkan dengan garis lurus (linear), maka regresi logistik akan sangat efektif dalam mengklasifikasikan data tersebut ke dalam dua kelompok berbeda. Dalam konteks ini, variabel target (y) hanya memiliki dua nilai, sehingga model bisa menentukan batas yang jelas di antara keduanya.

Memberikan wawasan yang bermakna Regresi logistik juga bisa menunjukkan seberapa besar pengaruh suatu variabel terhadap hasil akhir, melalui koefisiennya. Selain itu, model ini juga menunjukkan apakah hubungan antar variabel bersifat positif atau negatif.

Persamaan dan Asumsi dalam Regresi Logistik Model ini menggunakan fungsi sigmoid, yaitu fungsi matematika berbentuk kurva S yang mengubah nilai output ke dalam rentang antara 0 dan 1. Jika hasil dari fungsi sigmoid (yang merepresentasikan probabilitas) lebih tinggi dari ambang batas tertentu, maka model akan mengklasifikasikan data sebagai bagian dari suatu kelas. Sebaliknya, jika di bawah ambang tersebut, maka data dianggap tidak termasuk dalam kelas itu.

Selain itu, jika keluaran dari fungsi sigmoid (yaitu probabilitas yang diperkirakan) lebih besar dari ambang batas yang telah ditentukan dalam grafik, maka model akan memprediksi bahwa suatu observasi termasuk dalam kelas tertentu. Sebaliknya, jika nilai probabilitas tersebut lebih kecil dari ambang batas, model akan memprediksi bahwa observasi tersebut tidak termasuk ke dalam kelas tersebut

Sebagai contoh:

-Jika hasil fungsi sigmoid lebih dari 0,5, maka output dianggap sebagai 1 (kelas positif). -Jika hasilnya kurang dari 0,5, maka output diklasifikasikan sebagai 0 (kelas negatif). -Jika grafik menuju ke arah negatif secara ekstrem, maka nilai prediksi y akan menjadi 0,dan sebaliknya.

Fungsi sigmoid digunakan dalam regresi logistik untuk mengubah nilai prediksi menjadi probabilitas. Rumus dari fungsi sigmoid yaitu: \[ f(x) = \frac{1}{1 + e^{-x}} \]

# Simulasi data untuk regresi logistik 
set.seed(42)
n <- 100 
x <- seq(-4, 4, length.out = n)
log_odds <- -0.5 + 1.5 * x 
prob <- 1 / (1 + exp(-log_odds))
y <- rbinom(n, 1, prob)
# Buat data frame 
data <- data.frame(x=x, y=y, prob = prob)

Plot Kurva Sigmoid

# Visualisasi menggunakan base R plot
# Visualisasi menggunakan base R
# Plot the data
plot(x, y, 
     pch = 16, 
     col = "gray60",
     xlab = "X", 
     ylab = "Y / Probabilitas", 
     main = "Simulasi Regresi Logistik dengan Kurva Sigmoid")
lines(x, prob, col = "blue", lwd = 2)
abline(h = 0.5, col = "red", lty = 2)

legend("topleft", legend = c("Data Biner (0/1)", "Kurva Logistik", "Ambang 0.5"), 
       col = c("gray60", "blue", "red"), 
       pch = c(16, NA, NA), 
       lty = c(NA, 1, 2), 
       lwd = c(NA, 2, 1), 
       pt.cex = 1.5, bty = "n")

Kurva sigmoid dalam regresi logistik menunjukkan hubungan non-linear antaravariabel prediktor dan probabilitas output. Pendekatan ini efektif untuk klasifikasi biner seperti deteksi penyakit, kelulusan siswa, dan prediksi ya/tidak.

Fungsi Link (Logit): Fungsi link logit dapat dinyatakan dalam bentuk berikut:

\[ g(\mu) = \log \left( \frac{\mu}{1 - \mu} \right) \]

Model Regresi Logistik:

\[ \log \left( \frac{\mu}{1 - \mu} \right) = X\beta \]

Fungsi Inverse Link:

\[ \mu = \frac{\exp(X\beta)}{1 + \exp(X\beta)} \]

Estimasi Parameter

Metode estimasi parameter pada GLM umumnya menggunakan Maximum Likelihood Estimation (MLE).

Log-likelihood fungsi untuk regresi logistik:

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

\[ \pi_i = \frac{\exp(X\beta)}{1 + \exp(X\beta)} \]

Contoh Kasus dengan R

Misalkan kita memiliki:

#Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)

#Buat data biner
game <- rep(1:23, Attempts)
result <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))

#Gabung jadi data frame
data <- data.frame(Game = game, Result = result)

#Definisikan x dan y
x <- data$Game
y <- data$Result
head(data)

Estimasi Regresi Logistik Estimasi parameter model regresi logistik dapat menggunakan ‘glm’ function

model <- glm(y ~ x, data = data, family = binomial) 
summary(model)
## 
## Call:
## glm(formula = y ~ x, family = binomial, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept)  0.02085    0.26206   0.080    0.937
## x           -0.01546    0.01845  -0.838    0.402
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 408.06  on 295  degrees of freedom
## Residual deviance: 407.35  on 294  degrees of freedom
## AIC: 411.35
## 
## Number of Fisher Scoring iterations: 3

Interpretasi Koefisien

• Intercept ((Intercept)): nilai log-odds saat x=0

• Koefisien x: perubahan log-odds untuk setiap satu satuan peningkatan pada x

• Nilai p: menunjukkan signifikansi statistik dari prediktor

Koefisien log-odds dapat ditransformasikan ke odds ratio:

exp(coef(model))
## (Intercept)           x 
##   1.0210698   0.9846578

Prediksi dan Visualisasi Kurva Logit

library(ggplot2)
data$pred <- predict(model, type = "response")
ggplot(data, aes(x=x, y=y)) +
  geom_point(alpha = 0.5, color = "gray40") +
  geom_line(aes(y=pred), color = "blue", linewidth = 1.5) +
  labs(title = "Kurva Logit pada Regresi Logistik",
       x="X (Prediktor)",
       y="Probabilitas / Respons") + 
theme_minimal()

Evaluasi Model

Kita dapat melihat akurasi model berdasarkan klasifikasi threshold 0.5

# Klasifikasi dan confusion matrix 
data$pred_class <- ifelse(data$pred > 0.5, 1, 0) 
table(Predicted = data$pred_class, Actual = y)
##          Actual
## Predicted   0   1
##         0 160 131
##         1   1   4

GLM adalah kerangka model fleksibel untuk berbagai jenis data dan distribusi. Regresi logistik merupakan salah satu contoh penting dari GLM, sangat berguna dalam analisis data kategorik biner. Estimasi parameter dilakukan melalui metode MLE dan dapat diselesaikan secara efisien dengan fungsi glm di R.

Model Regresi Poisson

Regresi Poisson digunakan ketika variabel respons adalah data cacah (count data),yaitu bilangan bulat non negatif. Model ini merupakan bagian dari Generalized Linear Model(GLM) dengan asumsi bahwa distribusi variabel respons adalah distribusi Poisson.

Distribusi Poisson memiliki fungsi probabilitas:

\[ P(Y = y) = \frac{e^{-\lambda} \lambda^y}{y!} \]

Distribusi Poisson dapat dituliskan dalam bentuk exponential family sebagai berikut:

\[ f(y; \theta) = \exp \left\{ y \log(\lambda) - \lambda - \log(y!) \right\} \]

dengan:

  • \[ \theta = \log(\lambda) \]
  • \[ b(\theta) = e^{\theta} = \lambda \]
  • \[ \phi = 1 \]
  • \[ c(y, \phi) = -\log(y!) \]

Maka distribusi Poisson termasuk dalam exponential family.

Fungsi Link

Fungsi link kanonik untuk distribusi Poisson adalah fungsi logaritma:

\[ g(\mu) = \log(\mu) \]

Sehingga modelnya menjadi:

\[ \log(\mu_i) = x_i^T \beta \]

dan fungsi inverse link:

\[ \mu_i = \exp(x_i^T \beta) \]

Estimasi Parameter

Estimasi parameter \(( \beta )\) dilakukan dengan metode Maximum Likelihood Estimation (MLE).
Log-likelihood fungsi untuk regresi Poisson:

\[ l(\beta) = \sum_{i=1}^n \left[ y_i x_i^T \beta - \exp(x_i^T \beta) - \log(y_i!) \right] \]

Nilai \(( \beta )\) dapat diperoleh melalui metode numerik seperti iterasi Newton-Raphson.

Contoh

Misalkan kita memiliki data sebagai berikut:

# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)

# Buat data biner
game <- rep(1:23, Attempts)
result <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))
x <- data$Game
y <- data$Result

# Gabung jadi data frame
data <- data.frame(x, y)
head(data)
lambda <- exp(0.3 + 0.6 * x)

Estimasi Regresi Poisson

poisson_model <- glm(y ~ x, data = data, family = poisson())
summary(poisson_model)
## 
## Call:
## glm(formula = y ~ x, family = poisson(), data = data)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -0.679572   0.189073  -3.594 0.000325 ***
## x           -0.008375   0.013544  -0.618 0.536344    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 211.97  on 295  degrees of freedom
## Residual deviance: 211.59  on 294  degrees of freedom
## AIC: 485.59
## 
## Number of Fisher Scoring iterations: 5

Plot Hasil Prediksi

plot(x, y, pch = 16, col = "darkgray", main = "Data dan Hasil Prediksi") 
newdata <- data.frame(x=seq(min(x), max(x), length.out = 100)) 
pred <- predict(poisson_model, newdata = newdata, type = "response")
lines(newdata$x, pred, col = "blue", lwd = 2)

Diagnostik dan Overdispersion

Salah satu asumsi penting dari model Poisson adalah bahwa mean dan varians dari variabel respons adalah sama (𝐸[𝑌]=𝑉𝑎𝑟[𝑌]). Jika varians lebih besar dari mean, maka terjadi overdispersion. Untuk mendeteksi overdispersion:

dispersion <- sum(residuals(poisson_model, type = "pearson")^2) / poisson_model$df.residual 
dispersion
## [1] 0.54753

Jika nilai dispersion > 1, maka overdispersion mungkin terjadi dan model alternatif seperti Negative Binomial Regression dapat digunakan.

exp(coef(poisson_model))
## (Intercept)           x 
##   0.5068336   0.9916599

Visualisasi Prediksi

data$predicted <- predict(poisson_model, type = "response")
library(ggplot2) 
ggplot(data, aes(x=x, y=y)) +
  geom_jitter(width = 0.2, alpha = 0.6) + 
  geom_point(aes(y=predicted), shape = 18, size = 3, color = "black") +
  labs(title = "Prediksi", x="Game", y="Result") + 
  theme_minimal()

Evaluasi Model

# Plot residuals 
plot(poisson_model$residuals,
     main = "Residual Plot",
     ylab = "Residual",
     xlab = "Index",
     pch = 19,
     col = "steelblue")

abline(h = 0, col = "red", lty = 2)

Inferensi GLM

Dalam Generalized Linear Model (GLM), inferensi statistik membutuhkan pemahaman terhadap ekspektasi dan variansi dari estimator model, terutama untuk mengembangkan alat-alat uji seperti Wald test, Likelihood Ratio test, dan interval kepercayaan.

Ekspektasi dan Varians dalam GLM

  1. Ekspektasi Estimator

Ekspektasi menunjukkan apakah suatu estimator tak bias, yaitu:

\[ \mathbb{E}[\hat{\beta}] = \beta \]

Dalam GLM, MLE dari \(\hat{\beta}\) bersifat asymptotically unbiased.

  1. Varians Estimator

Varians menunjukkan presisi dari estimasi parameter:

\[ \text{Var}(\hat{\beta}) \approx (\mathbf{X}^T \mathbf{W} \mathbf{X})^{-1} \]

dimana \(\mathbf{W}\) adalah matriks bobot yang tergantung pada distribusi dan fungsi link.

Distribusi Asimptotik Estimator

Dengan ukuran sampel besar:

\[ \hat{\beta} \sim \mathcal{N}(\beta, \text{Var}(\hat{\beta})) \]

Distribusi ini adalah dasar dari:

  • Uji Wald
  • Confidence interval
  • P-value

Varians dalam GLM Tidak Konstan

Tidak seperti regresi linear (OLS) yang mengasumsikan homoskedastisitas:

\[ \text{Var}(Y_i) = \sigma^2 \]

Dalam GLM:

\[ \text{Var}(Y_i) = \phi V(\mu_i) \]

  • \(\phi\) = parameter dispersi
  • \(V(\mu)\) = fungsi varians

Contoh:

  • Poisson: \(V(\mu) = \mu\)
  • Binomial: \(V(\mu) = \mu(1 - \mu)\)

Contoh Regresi Poisson

# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)

# Buat data biner
x <- rep(1:23, Attempts)
y <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))

# Gabung jadi data frame
data <- data.frame(x = x, y = y)

# Ringkasan model regresi Poisson
model <- glm(y ~ x, data = data, family = poisson())
summary(model)
## 
## Call:
## glm(formula = y ~ x, family = poisson(), data = data)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -0.679572   0.189073  -3.594 0.000325 ***
## x           -0.008375   0.013544  -0.618 0.536344    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 211.97  on 295  degrees of freedom
## Residual deviance: 211.59  on 294  degrees of freedom
## AIC: 485.59
## 
## Number of Fisher Scoring iterations: 5
  • Ekspektasi: \(( \mathbb{E}[Y_i] = \mu_i )\)
  • Varians: \(( \text{Var}(Y_i) = \mu_i )\)

Kesimpulan:

  • Ekspektasi digunakan untuk mengetahui ketakbiasan estimasi
  • Varians digunakan untuk mengukur presisi dan menyusun uji statistik
  • Distribusi asimptotik dari \(( \hat{\beta} )\) sangat bergantung pada kedua konsep ini
  • Dalam GLM, varians sangat tergantung pada bentuk distribusi eksponensial dari data

Mencari Ekspektasi dan Varians dalam GLM

Ekspektasi

Jika diturunkan berdasarkan fungsi momen:

\[ E(Y) = \int y f(y; \theta) \, dy = \mu \]

Untuk keluarga eksponensial:

\[ \log f(y; \theta) = a(y) + b(\theta) y + c(\theta) \]

atau:

\[ \log f(y; \theta) = y \theta - b(\theta) + c(y) \]

Maka ekspektasi turunan pertama:

\[ U(\theta) = \frac{\partial \ell}{\partial \theta} = y - b'(\theta) \] Dan ekspektasi turunan pertama:

\[ \mathbb{E}[U(\theta)] = \mathbb{E}[y - b'(\theta)] = \mu - b'(\theta) = 0 \]

Maka:

\[ \mu = b'(\theta) \]

Varians
Turunan kedua:

\[ \frac{\partial^2 \ell}{\partial \theta^2} = -b''(\theta) \]

Sehingga:

\[ \text{Var}(Y) = b''(\theta) = \phi V(\mu) \]

Metode Penaksiran Parameter

Maximum Likelihood Estimation (MLE)

  • Prinsip dasar: memaksimumkan fungsi likelihood/log-likelihood.
  • Langkah:
    • Turunan pertama = 0
    • Turunan kedua < 0

Namun, karena bentuk GLM tidak eksplisit, digunakan metode numerik.

Metode Optimisasi: Newton-Raphson

  • Menggunakan score vector (gradien)
  • Menggunakan Hessian matrix

Iterasi:

\[ \beta^{(t+1)} = \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \] Fisher Scoring

  • Modifikasi Newton-Raphson, mengganti Hessian dengan matriks informasi Fisher.

IRLS (Iteratively Reweighted Least Square) - Modifikasi dari Fisher scoring, hasil estimasi mirip dengan Least Square.

Implementasi Newton-Raphson Statistik score ke-𝑗 Statistik score ke-𝑗:

\[ U_j(\beta) = \frac{\partial \log L(\beta)}{\partial \beta_j} \] Turunan kedua:

\[ H_{jk}(\beta) = \frac{\partial^2 \log L(\beta)}{\partial \beta_j \partial \beta_k} \]

Taylor expansion:

\[ U(\beta^*) \approx U(\beta) + H(\beta)(\beta^* - \beta) \]

Estimasi parameter:

\[ \hat{\beta} \approx \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \]

Diagnostik Model GLM

Diagnostik digunakan untuk mengevaluasi apakah model sudah tepat.

  • Uji formal
  • Grafik antara nilai prediksi vs nilai aktual

Statistik Devians

  • Mengukur apakah ada model lain yang lebih baik.
  • Nilai devians besar → model tidak cocok.
  • Devians adalah:

\[ D = 2 \sum \left[ y_i \log \left( \frac{y_i}{\hat{\mu}_i} \right) - (y_i - \hat{\mu}_i) \right] \]

  • Devians membandingkan model terhadap saturated model.
  • Devians kecil → model lebih cocok.

Statistik Chi-Kuadrat Pearson

  • Menguji apakah model lebih baik daripada tidak ada model sama sekali.
  • Statistik:

\[ X^2 = \sum \left( \frac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i} \right) \] • Jika signifikan → model lebih baik daripada tanpa model.

Catatan

• Untuk data yang dikelompokkan, statistik devians dan chi-kuadrat Pearson mengikuti distribusi ChiSquare.

• Untuk data tidak dikelompokkan, tidak mengikuti distribusi Chi-Square.

• Devians diminimalkan oleh MLE → cocok digunakan untuk evaluasi model.

Analisis Residual

• Residual adalah selisih antara observasi dengan prediksi.

• Dapat digunakan untuk memeriksa penyimpangan sistematis.

• Dapat diplot untuk menilai asumsi model

Detail Metode Estimasi dan Inferensi Regresi Logistik

Regresi logistik digunakan untuk memodelkan probabilitas dari variabel respons biner (0/1) berdasarkan satu atau lebih variabel prediktor. Estimasi parameter dilakukan menggunakan Maximum Likelihood Estimation (MLE) karena model tidak linear dalam parameternya.

Fungsi model logistik:

\[ \pi(x) = \frac{\exp(\beta_0 + \beta_1 x)}{1 + \exp(\beta_0 + \beta_1 x)} \]

Log-likelihood untuk \(( n )\) observasi:

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

Metode Newton-Raphson digunakan untuk mencari nilai parameter \(( \beta )\) yang memaksimalkan fungsi log-likelihood pada model regresi logistik.

Fungsi Log-Likelihood

Model regresi logistik untuk probabilitas:

\[ \pi_i = \frac{1}{1 + \exp(-x_i^\top \beta)} \]

Log-likelihood untuk \(( n )\) observasi:

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

Langkah-langkah Newton-Raphson Turunan dalam Estimasi Logistik

  1. Turunan Pertama (Score Function)

\[ U(\beta) = \frac{\partial \ell(\beta)}{\partial \beta} = X^\top (y - \pi) \]

  1. Turunan Kedua (Hessian Matrix)

\[ H(\beta) = -X^\top W X, \quad \text{dengan } W = \text{diag}(\pi_i (1 - \pi_i)) \]

  1. Iterasi Newton-Raphson

\[ \beta^{(t+1)} = \beta^{(t)} + (X^\top W^{(t)} X)^{-1} X^\top (y - \pi^{(t)}) \]

Estimasi MLE dengan Newton-Raphson (Manual di R)

# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)

# Buat data biner
x <- rep(1:23, Attempts)
y <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))

# Gabung jadi data frame
data <- data.frame(x = x, y = y)

beta_true <- c(-1, 2)
X <- cbind(1, x)
eta <- X %*% beta_true
p <- 1 / (1 + exp(-eta))

Newton-Raphson Iterasi Manual

beta <- c(0, 0)
tol <- 1e-6 
max_iter <- 100 
for (i in 1:max_iter) {
  eta <- X %*% beta
  pi_hat <- 1 / (1 + exp(-eta)) 
  W <- diag(as.numeric(pi_hat * (1 - pi_hat)))
  z <- eta + solve(W) %*% (y - pi_hat) 
  beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z 
  if (sum(abs(beta_new - beta)) < tol) break
  beta <- beta_new }
beta # hasil estimasi akhir
##          [,1]
##    0.02085078
## x -0.01546112

• Estimasi parameter pada model regresi logistik dilakukan dengan MLE. • Newton-Raphson adalah metode numerik yang digunakan untuk memaksimalkan log-likelihood.

• Iterasi didasarkan pada turunan pertama (score) dan kedua (Hessian).

• Prosedur identik dengan IRLS (Iteratively Reweighted Least Squares) dalam implementasi GLM.

Inferensi Parameter

  1. Uji Wald

Tujuan Uji Wald
Untuk menguji signifikansi parameter \(( \beta_j )\) dalam model regresi logistik:

  • \(( H_0: \beta_j = 0 )\)
  • \(( H_1: \beta_j \neq 0 )\)

Teori Wald Test
Dari teori estimasi MLE, estimator \(( \hat{\beta}_j )\) mendekati distribusi normal:

\[ \hat{\beta}_j \sim N(\beta_j, \text{Var}(\hat{\beta}_j)) \]

Jika \(( H_0 )\) benar (yaitu \(( \beta_j = 0 )\)), maka:

\[ Z = \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)} \sim N(0,1) \]

Dengan Statistik Wald:

\[ W = Z^2 = \left( \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)} \right)^2 \sim \chi^2_1 \]

Uji Wald Langkah demi Langkah

# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)

# Buat data biner
x <- rep(1:23, Attempts)
y <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))

# Gabung jadi data frame
data <- data.frame(x = x, y = y)

log_odds <- -0.5 + 1.2 * x
p <- 1 / (1 + exp(-log_odds))

model <- glm(y ~ x, data = data, family = binomial)

summary(model)
## 
## Call:
## glm(formula = y ~ x, family = binomial, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept)  0.02085    0.26206   0.080    0.937
## x           -0.01546    0.01845  -0.838    0.402
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 408.06  on 295  degrees of freedom
## Residual deviance: 407.35  on 294  degrees of freedom
## AIC: 411.35
## 
## Number of Fisher Scoring iterations: 3

Langkah 1: Ambil nilai koefisien dan standard error (SE)

beta_hat <- coef(model)["x"]
se_beta <- summary(model)$coefficients["x", "Std. Error"]

Langkah 2: Hitung Statistik Z

Z <- beta_hat / se_beta
Z
##          x 
## -0.8380593

Langkah 3: Hitung Statistik Wald

Wald_stat <- Z^2
Wald_stat
##         x 
## 0.7023434

Langkah 4: Hitung p-value

p_value <- 1 - pchisq(Wald_stat, df = 1)
p_value
##         x 
## 0.4019974

Interpretasi

• Jika p-value < 0.05, maka koefisien signifikan → variabel prediktor berpengaruh.

• Jika p-value > 0.05, maka tidak ada cukup bukti untuk menolak - \(( H_0 )\)

Kesimpulan

Uji Wald didasarkan pada rasio antara estimasi parameter dan standar error-nya. Dengan menaikkan nilai Z menjadi kuadrat (Z²), kita memperoleh distribusi Chi-Square untuk pengujian hipotesis parameter individual dalam model regresi logistik

  1. Uji Likelihood Ratio (Chi-Square)

Bandingkan model penuh dengan model tanpa prediktor

# Model null 
model_null <- glm(y ~ 1, data = data, family = binomial)

# Likelihood ratio test 
anova(model_null, model, test = "Chisq")

Evaluasi Kebaikan Model

  1. Akaike Information Criterion (AIC)

Semakin kecil AIC, semakin baik model.

AIC(model)
## [1] 411.3528
  1. Bayesian Information Criterion (BIC)

Alternatif terhadap AIC, menghukum kompleksitas model.

BIC(model)
## [1] 418.7335

• Estimasi regresi logistik dilakukan dengan MLE melalui iterasi Newton-Raphson.

• Inferensi parameter dapat dilakukan dengan uji Wald dan likelihood ratio (uji Chi-Square).

• AIC dan BIC digunakan untuk mengevaluasi kompleksitas dan kecocokan model.

Detail Metode Estimasi dan Inferensi Regresi Poisson

Model regresi Poisson digunakan untuk memodelkan data count (jumlah kejadian) dimana variabel respons mengikuti distribusi Poisson. Estimasi dilakukan dengan Maximum Likelihood Estimation (MLE), dan inferensi dilakukan dengan uji Wald dan Likelihood Ratio Test.

Model Regresi Poisson

Distribusi Poisson:

\[ P(Y_i = y_i) = \frac{e^{-\lambda_i} \lambda_i^{y_i}}{y_i!} \]

Model regresi Poisson:

\[ \log(\lambda_i) = x_i^T \beta \] Estimasi Parameter (MLE)

Log-likelihood fungsi:

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

Dengan:

\[ \lambda_i = \exp(x_i^T \beta) \]

Estimasi dilakukan dengan metode iterasi (IRLS)

Estimasi parameter model regresi Poisson menggunakan pendekatan
Maximum Likelihood Estimation (MLE) dengan metode
Iteratively Reweighted Least Squares (IRLS) secara manual, tanpa menggunakan glm().

Tahap 1: Definisikan Model Regresi Poisson

\[ \log(\lambda_i) = x_i^T \beta \quad \Rightarrow \quad \lambda_i = \exp(x_i^T \beta) \]

Tahap 2: Mencari Log-likelihood yang dimaksimumkan

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

\[ \beta^{(t+1)} = \left( \mathbf{X}^T \mathbf{W}^{(t)} \mathbf{X} \right)^{-1} \mathbf{X}^T \mathbf{W}^{(t)} \mathbf{z}^{(t)} \]

Dengan:

  • \[\mathbf{W} = \text{diag}(\lambda_i)\]
  • \[\mathbf{z} = \eta + \frac{y - \lambda}{\lambda}\]

dan

\[ \eta_i = \log(\lambda_i) = x_i^T \beta \] Simulasi Data

# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)

# Buat data biner
x <- rep(1:23, Attempts)
y <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))

# Gabung jadi data frame
data <- data.frame(x = x, y = y)

X <- cbind(1, x)  # Tambah intercept
beta_true <- c(0.5, 0.8)
eta <- X %*% beta_true
lambda <- exp(eta)

IRLS Manual Step-by-Step

# Inisialisasi 
beta <- c(0, 0) 
tol <- 1e-6 
max_iter <- 100 

for (i in 1:max_iter) {
  eta <- X %*% beta 
  lambda <- exp(eta) 
  W <- diag(as.numeric(lambda))
  z <- eta + (y - lambda) / lambda
  
  beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z 
  
  if (sum(abs(beta_new - beta)) < tol) { 
    cat("Konvergen pada iterasi ke-", i, "\n")
    break 
    } 
  beta <- beta_new
}
## Konvergen pada iterasi ke- 5
beta # hasil estimasi
##           [,1]
##   -0.679572484
## x -0.008375021

Perbandingan dengan glm

model_glm <- glm(y ~ x, family = poisson()) 
summary(model_glm)
## 
## Call:
## glm(formula = y ~ x, family = poisson())
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -0.679572   0.189073  -3.594 0.000325 ***
## x           -0.008375   0.013544  -0.618 0.536344    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 211.97  on 295  degrees of freedom
## Residual deviance: 211.59  on 294  degrees of freedom
## AIC: 485.59
## 
## Number of Fisher Scoring iterations: 5

• IRLS memberikan cara iteratif untuk menghitung estimasi MLE dalam regresi Poisson.

• Hasil manual IRLS sangat mendekati hasil glm() dari R.

• Metode ini memberikan pemahaman mendalam atas mekanisme di balik fungsi glm().

Pengujian hipotesis Uji Wald

Untuk menguji H0: = 0

# Koefisien dan standar error
coef_val <- coef(model)[2] 
se_val <- summary(model)$coefficients[2, 2]
wald_z <- coef_val / se_val 
wald_chisq <- wald_z^2 
p_value <- 1 - pchisq(wald_chisq, df = 1)
cat("Z:", wald_z, "\nChi-Square:", wald_chisq, "\np-value:", p_value)
## Z: -0.8380593 
## Chi-Square: 0.7023434 
## p-value: 0.4019974

Uji Likelihood Ratio (Chi-Square)

model_null <- glm(y ~ 1, family = poisson(), data = data)
anova(model_null, model, test = "Chisq")

Evaluasi Model (AIC & BIC)

AIC(model)
## [1] 411.3528
BIC(model)
## [1] 418.7335

• Estimasi parameter regresi Poisson dilakukan menggunakan MLE

• Uji Wald dan Likelihood Ratio digunakan untuk pengujian hipotesis

• AIC dan BIC digunakan untuk evaluasi dan pemilihan model terbaik

Regresi Logistik dengan Prediktor Nominal, Ordinal, dan Rasio

Regresi logistik merupakan salah satu metode analisis statistik yang digunakan untuk memodelkan hubungan antara satu variabel respons biner (dua kategori) dengan satu atau lebih variabel prediktor. Dalam penerapannya, prediktor yang digunakan dapat memiliki skala pengukuran berbeda, yaitu:

Nominal: Variabel yang tidak memiliki urutan logis antar kategorinya, hanya sebagai pembeda.

Contoh: jenis kelamin (laki-laki, perempuan). Dalam analisis regresi logistik, data nominal diubah menjadi variabel dummy (misal: Female = 1, Male = 0).

Ordinal: Variabel yang memiliki urutan logis antar kategori, tetapi jarak antar kategori belum tentu sama.

Contoh: tingkat pendidikan (SMA < Sarjana < Master < Doktor). Dalam analisis, variabel ordinal bisa diperlakukan:

– Sebagai nominal dengan membuat dummy untuk setiap kategori.

– Seperti rasio dengan memberi nilai tertentu (misalnya SMA = 1, Sarjana = 2, Master = 3, Doktor = 4) untuk mencerminkan urutan dan jarak antar tingkatannya.

Rasio: Variabel numerik kontinu yang memiliki nol absolut dan rasio bermakna.

Contoh: pendapatan (dalam juta rupiah per bulan). Data ini dapat langsung digunakan dalam model tanpa transformasi khusus.

Simulasi Data

Dalam dunia bisnis modern, startup merupakan salah satu bentuk usaha yang paling dinamis namun juga paling rentan mengalami kegagalan. Untuk memahami faktor-faktor yang memengaruhi keberhasilan atau kegagalan sebuah startup, dilakukan simulasi data terhadap 250 proyek startup dari berbagai bidang. Setiap proyek dikategorikan berdasarkan bidang usaha (field), tingkat pengalaman founder (experience), serta jumlah modal awal (capital). Tujuan analisis ini adalah untuk memodelkan probabilitas kegagalan proyek menggunakan regresi logistik biner, di mana kegagalan didefinisikan sebagai 1 dan keberhasilan sebagai 0.

Variabel-variabel yang digunakan dalam studi ini adalah:

Failure (biner): Status hasil proyek startup. Nilai 1 menunjukkan startup gagal, sedangkan nilai 0 menunjukkan startup berhasil.

Field (kategori): Bidang usaha startup, terdiri atas empat kategori, yaitu Fintech, Healthtech, Edtech, dan E-commerce.

Experience (kategori): Tingkat pengalaman pendiri (founder) startup, dikategorikan menjadi Low, Medium, dan High.

Capital (numerik): Modal awal startup yang dimiliki pada saat pendirian, dinyatakan dalam ribuan dolar Amerika Serikat (USD).

Dengan menggunakan ketiga variabel prediktor tersebut, probabilitas kegagalan setiap proyek dianalisis untuk melihat bagaimana bidang usaha, pengalaman pendiri, dan besar modal awal berkontribusi terhadap risiko kegagalan startup

Simulasi dataset akan dibuat menggunakan 250 observasi

library(tibble)
library(dplyr)
set.seed(123)

# Simulasi variabel
n <- 250

# Variabel X: bidang startup
field <- sample(c("Fintech", "Healthtech", "Edtech", "E-commerce"), 
                n, replace = TRUE, prob = c(0.3, 0.2, 0.2, 0.3))

# Variabel Y: tingkat pengalaman founder
experience <- sample(c("Low", "Medium", "High"), 
                     n, replace = TRUE, prob = c(0.4, 0.4, 0.2))

# Variabel tambahan: jumlah modal awal (dalam ribu USD)
capital <- rnorm(n, mean = 500, sd = 100)

# Membuat peluang kegagalan
logit_p <- -1.5 + 
           0.7 * (field == "E-commerce") + 
           0.5 * (experience == "Low") + 
           0.002 * (500 - capital)  # semakin kecil modal, semakin besar peluang gagal

p <- 1 / (1 + exp(-logit_p))

set.seed(123)
failure <- rbinom(n, 1, p)  # 1 = gagal, 0 = berhasil

# Gabungkan menjadi data frame
library(tibble)
startup_data <- tibble(failure, field, experience, capital)

head(startup_data)

Eksplorasi Data

startup_data %>%
  dplyr::group_by(failure) %>%
  dplyr::summarise(
    Jumlah = dplyr::n(),
Rata2_Modal = mean(capital)
)

Perlakuan Variabel Ordinal

Treat Sebagai Nominal (Dummy)

startup_data_nominal <- startup_data %>%
  mutate(
    field = factor(field, levels = c("Fintech", "Healthtech", "Edtech", "E-commerce")),
    experience = factor(experience, levels = c("Low", "Medium", "High"))
  )

model_nominal <- glm(failure ~ field + experience + capital,
                     data = startup_data_nominal,
                     family = binomial)
summary(model_nominal)
## 
## Call:
## glm(formula = failure ~ field + experience + capital, family = binomial, 
##     data = startup_data_nominal)
## 
## Coefficients:
##                    Estimate Std. Error z value Pr(>|z|)   
## (Intercept)      -1.260e+01  1.848e+03  -0.007  0.99456   
## fieldHealthtech   2.659e+01  1.848e+03   0.014  0.98852   
## fieldEdtech       1.824e+01  1.848e+03   0.010  0.99212   
## fieldE-commerce   1.758e+01  1.848e+03   0.010  0.99241   
## experienceMedium -3.280e+00  1.157e+00  -2.835  0.00458 **
## experienceHigh   -5.289e+00  1.703e+00  -3.107  0.00189 **
## capital          -1.422e-02  4.432e-03  -3.208  0.00134 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 260.899  on 249  degrees of freedom
## Residual deviance:  59.771  on 243  degrees of freedom
## AIC: 73.771
## 
## Number of Fisher Scoring iterations: 19

Interpretasi Koefisien

Intercept (-12.60)

-Ini adalah log-odds dasar untuk startup yang bergerak di bidang Fintech, dengan founder berpengalaman Low, dan modal standar (capital = 0).

-Nilai p = 0.99456 (tidak signifikan di level 5%), artinya baseline ini tidak berbeda signifikan terhadap probabilitas gagal.

-Odds Ratio = exp(-12.60) ≈ 0.00000336 → Sangat kecil. Peluang kegagalan pada kondisi baseline hampir nol, meski tidak signifikan.

fieldHealthtech (+26.59)

-Startup di bidang Healthtech memiliki log-odds kegagalan lebih tinggi sebesar 26.59 dibandingkan Fintech (referensi).

-p = 0.98852 (tidak signifikan), artinya tidak ada bukti signifikan bahwa Healthtech lebih berisiko gagal dibanding Fintech.

-Odds Ratio = exp(26.59) ≈ 3.46 × 10¹¹ → Secara matematis sangat besar, tapi tidak signifikan.

fieldEdtech (+18.24)

-Startup di bidang Edtech memiliki log-odds kegagalan lebih tinggi sebesar 18.24 dibandingkan Fintech.

-p = 0.99212 (tidak signifikan), sehingga perbedaan tidak cukup bukti untuk menyatakan risiko kegagalan Edtech lebih tinggi dibanding Fintech.

-Odds Ratio = exp(18.24) ≈ 9.2 × 10⁷ → Besar, tapi tidak signifikan.

fieldE-commerce (+17.58)

-Startup di bidang E-commerce memiliki log-odds kegagalan lebih tinggi sebesar 17.58 dibandingkan Fintech.

-p = 0.99241 (tidak signifikan), sehingga tidak cukup bukti bahwa E-commerce lebih rentan gagal dibanding Fintech.

-Odds Ratio = exp(17.58) ≈ 4.3 × 10⁷ → Besar, namun tidak signifikan.

experienceMedium (-3.28)

-Founder dengan pengalaman Medium memiliki log-odds kegagalan lebih rendah sebesar 3.28 dibandingkan Low (referensi).

-p = 0.00458 (signifikan di level 5%), berarti pengalaman Medium signifikan menurunkan peluang kegagalan.

-Odds Ratio = exp(-3.28) ≈ 0.038 → Peluang gagal hanya sekitar 3.8% dibandingkan dengan founder berpengalaman Low.

experienceHigh (-5.289)

-Founder dengan pengalaman High memiliki log-odds kegagalan lebih rendah sebesar 5.289 dibandingkan Low.

-p = 0.00189 (signifikan di level 5%), menunjukkan bahwa pengalaman tinggi sangat signifikan mengurangi kegagalan.

-Odds Ratio = exp(-5.289) ≈ 0.005 → Peluang gagal hanya sekitar 0.5% dari peluang founder berpengalaman rendah.

capital (-0.01422)

-Setiap kenaikan modal sebesar 1 unit (asumsinya dalam satuan modal standar), menurunkan log-odds kegagalan sebesar 0.01422.

-p = 0.00134 (signifikan di level 5%), artinya modal lebih besar signifikan menurunkan peluang kegagalan.

-Odds Ratio = exp(-0.01422) ≈ 0.9859 → Setiap tambahan modal 1 unit mengurangi peluang gagal sekitar 1.4%.

Interpretasi Goodness-of-Fit

Null Deviance (260.899): Deviance model tanpa prediktor.

Residual Deviance (59.771): Deviance model dengan prediktor. → Penurunan dari null deviance ke residual deviance menunjukkan bahwa model membawa informasi yang cukup baik, meski penurunannya sedang saja (tidak drastis).

AIC (73.771): → Digunakan untuk membandingkan model: semakin kecil AIC, semakin baik model dalam menyeimbangkan goodness-of-fit dan kompleksitas. → Nilai AIC ini tergolong cukup rendah untuk data simulasi seperti ini.

Signifikansi Model

Variabel yang signifikan: → experienceMedium, experienceHigh, dan capital berpengaruh signifikan terhadap penurunan peluang gagal startup (p < 0.01).

Variabel yang tidak signifikan: → Variabel field (baik Healthtech, Edtech, maupun E-commerce) tidak signifikan mempengaruhi peluang kegagalan dibandingkan Fintech.

Arah hubungan: → Semakin tinggi pengalaman dan modal, semakin kecil peluang startup untuk mengalami kegagalan.

Kesimpulan Praktis

Pengalaman founder (medium atau high) merupakan faktor kunci dalam menurunkan risiko kegagalan startup.

Modal yang lebih besar memberikan perlindungan signifikan terhadap kemungkinan gagal.

Bidang usaha (field) startup tidak menunjukkan pengaruh yang signifikan terhadap kegagalan dalam model ini, menunjukkan bahwa faktor internal (seperti pengalaman dan modal) lebih berperan.

Model sudah lebih baik dibandingkan model null, tapi masih ada ruang untuk peningkatan model, misalnya dengan menambahkan variabel baru atau memperbaiki struktur data.

Treat Sebagai Rasio (Numeric Rank)

startup_data_numeric <- startup_data %>%
  mutate(
    experience_numeric = case_when(
      experience == "Low" ~ 1,
      experience == "Medium" ~ 2,
      experience == "High" ~ 3
    )
  )

model_numeric <- glm(failure ~ field + experience_numeric + capital, 
                     data = startup_data_numeric, 
                     family = binomial)

summary(model_numeric)
## 
## Call:
## glm(formula = failure ~ field + experience_numeric + capital, 
##     family = binomial, data = startup_data_numeric)
## 
## Coefficients:
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         7.684e+00  2.372e+00   3.240 0.001197 ** 
## fieldEdtech         6.730e-01  8.064e-01   0.835 0.403978    
## fieldFintech       -1.758e+01  1.861e+03  -0.009 0.992466    
## fieldHealthtech     9.160e+00  1.753e+00   5.224 1.75e-07 ***
## experience_numeric -2.799e+00  8.325e-01  -3.362 0.000773 ***
## capital            -1.417e-02  4.363e-03  -3.247 0.001168 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 260.899  on 249  degrees of freedom
## Residual deviance:  60.258  on 244  degrees of freedom
## AIC: 72.258
## 
## Number of Fisher Scoring iterations: 19

Interpretasi Koefisien

Intercept (+7.684)

Ini adalah log-odds dasar untuk startup pada bidang E-commerce (sebagai baseline), dengan pengalaman dan modal nol. Nilai ini signifikan (p = 0.0012), menunjukkan bahwa probabilitas kegagalan dasar cukup tinggi.

fieldEdtech (+0.6730)

Startup pada bidang Edtech memiliki log-odds kegagalan yang lebih tinggi sebesar 0.673 dibandingkan E-commerce, namun tidak signifikan (p = 0.4040). Ini menunjukkan bahwa perbedaan antara Edtech dan E-commerce dalam hal kegagalan tidak cukup kuat untuk disimpulkan secara statistik.

fieldFintech (–17.5800)

Startup pada bidang Fintech memiliki log-odds kegagalan yang jauh lebih rendah dibandingkan E-commerce, namun tidak signifikan (p = 0.9925). Nilai koefisien ekstrem ini kemungkinan muncul akibat data yang sangat tidak seimbang (jumlah Fintech sangat sedikit atau tidak gagal).

fieldHealthtech (+9.1600)

Startup pada bidang Healthtech memiliki log-odds kegagalan lebih tinggi secara signifikan sebesar 9.16 dibandingkan E-commerce. p < 0.001, sangat signifikan.

Odds Ratio ≈ exp(9.16) ≈ 9527.55 Artinya, peluang gagal startup Healthtech sekitar 9527 kali lebih tinggi dari E-commerce (namun interpretasi ini harus hati-hati karena nilainya sangat ekstrem dan bisa dipengaruhi data outlier/struktur data tertentu).

experience_numeric (–2.799)

Setiap kenaikan 1 tingkat pengalaman (misalnya dari Low → Medium, atau Medium → High) menurunkan log-odds kegagalan sebesar 2.799.

p = 0.0008 (sangat signifikan)

Odds Ratio = exp(–2.799) ≈ 0.061 Artinya, setiap peningkatan pengalaman menurunkan peluang gagal sekitar 94%.

capital (–0.01417)

Setiap penambahan 1 unit modal (misalnya 1 juta), menurunkan log-odds kegagalan sebesar 0.0142.

p = 0.0012 (sangat signifikan)

Odds Ratio = exp(–0.01417) ≈ 0.986 Artinya, setiap tambahan 1 juta modal menurunkan peluang gagal sekitar 1.4%.

Interpretasi Goodness-of-Fit

Null deviance (260.899): Deviance model tanpa prediktor.

Residual deviance (60.258): Deviance model dengan prediktor. → Penurunan dari null deviance ke residual deviance menunjukkan bahwa prediktor memberikan tambahan informasi yang sangat berarti.

AIC (72.258) → Digunakan untuk membandingkan model: semakin kecil AIC, semakin baik keseimbangan antara goodness-of-fit dan kompleksitas model.

Signifikansi Model

Variabel experience_numeric (p = 0.000773) dan capital (p = 0.001168) sangat signifikan dalam menurunkan peluang kegagalan.

Variabel fieldHealthtech (p = 1.75e-07) juga signifikan, dengan log-odds kegagalan jauh lebih tinggi dibandingkan bidang referensi (E-commerce).

Variabel fieldEdtech dan fieldFintech tidak signifikan secara statistik, artinya tidak terbukti berbeda secara nyata dari bidang E-commerce.

Kesimpulan Praktis

Semakin tinggi tingkat pengalaman, semakin kecil peluang startup mengalami kegagalan.

Modal (capital) yang lebih besar secara konsisten menurunkan peluang kegagalan, walaupun efeknya kecil (koefisien -0.014).

Startup di bidang Healthtech memiliki risiko kegagalan jauh lebih tinggi dibandingkan bidang E-commerce (log-odds = 9.16).

Variabel bidang lain seperti Edtech dan Fintech tidak menunjukkan perbedaan signifikan dalam hal peluang kegagalan.

Model ini menunjukkan performa yang sangat baik dalam menjelaskan kegagalan startup, dan masih bisa ditingkatkan dengan menambah variabel lain jika tersedia

Perbandingan Model

list(
AIC_Nominal = AIC(model_nominal),
AIC_Numeric = AIC(model_numeric)
)
## $AIC_Nominal
## [1] 73.77055
## 
## $AIC_Numeric
## [1] 72.25842

Keterangan: Semakin kecil nilai AIC, semakin bagus untuk dipakai

Goodness-of-Fit Model

nullmod <- glm(failure ~ 1, data = startup_data, family = binomial)
r2_nominal <- 1 - (logLik(model_nominal)/logLik(nullmod))
r2_numeric <- 1 - (logLik(model_numeric)/logLik(nullmod))
list(
McFadden_R2_Nominal = r2_nominal,
McFadden_R2_Numeric = r2_numeric
)
## $McFadden_R2_Nominal
## 'log Lik.' 0.7709056 (df=7)
## 
## $McFadden_R2_Numeric
## 'log Lik.' 0.7690357 (df=6)

Interpretasi: Semakin besar nilainya, semakin baik model memprediksi dibandingkan model null

Visualisasi Prediksi

library(ggplot2)
startup_data_nominal <- startup_data_nominal %>% 
  mutate(predicted = predict(model_nominal, type = "response"))

startup_data_numeric <- startup_data_numeric %>% 
  mutate(predicted = predict(model_numeric, type = "response"))

# Plot untuk model nominal
startup_data_nominal %>%
  ggplot(aes(x = capital, y = predicted, color = experience)) +
  geom_point(alpha = 0.6) +
  labs(
    title = "Prediksi Probabilitas (Ordinal sebagai Nominal)",
    x = "Capital",
    y = "Prediksi Probabilitas"
  ) +
  theme_minimal()

# Plot untuk model numeric
startup_data_numeric %>%
  ggplot(aes(x = capital, y = predicted, color = experience)) +
  geom_point(alpha = 0.6) +
  labs(
    title = "Prediksi Probabilitas (Ordinal sebagai Numeric)",
    x = "Capital",
    y = "Prediksi Probabilitas"
  ) +
  theme_minimal()

Ringkasan Koefisien Model

Ringkasan Koefisien Model Nominal (tanpa Odds Ratio)
term estimate std.error statistic p.value conf.low conf.high
(Intercept) -12.600 1847.902 -0.007 0.995 NA 118.728
fieldHealthtech 26.591 1847.901 0.014 0.989 -38.652 710.624
fieldEdtech 18.240 1847.901 0.010 0.992 -113.088 NA
fieldE-commerce 17.581 1847.901 0.010 0.992 -113.747 NA
experienceMedium -3.280 1.157 -2.835 0.005 -6.323 -1.395
experienceHigh -5.289 1.703 -3.107 0.002 -9.349 -2.472
capital -0.014 0.004 -3.208 0.001 -0.024 -0.006
Ringkasan Koefisien Model Nominal (tanpa Odds Ratio)
term estimate std.error statistic p.value conf.low conf.high
(Intercept) 7.684 2.372 3.240 0.001 3.441 12.907
fieldEdtech 0.673 0.806 0.835 0.404 -0.949 2.290
fieldFintech -17.575 1861.204 -0.009 0.992 NA 114.728
fieldHealthtech 9.160 1.753 5.224 0.000 6.323 13.388
experience_numeric -2.799 0.833 -3.362 0.001 -4.803 -1.425
capital -0.014 0.004 -3.247 0.001 -0.024 -0.006

Kesimpulan

-Semakin tinggi tingkat pengalaman, semakin kecil peluang startup mengalami kegagalan.

-Modal yang lebih besar memberikan perlindungan signifikan terhadap kemungkinan gagal.

-Bidang usaha (field) startup tidak menunjukkan pengaruh yang signifikan terhadap kegagalan dalam model ini, menunjukkan bahwa faktor internal (seperti pengalaman dan modal) lebih berperan.

Pemilihan Model Regresi Logistik dan Evaluasi

Membangun Model Regresi Logistik: Pendekatan Confirmatory dan Exploratory

Dalam analisis regresi logistik, pemilihan model sangat krusial untuk mendapatkan model yang baik dalam memprediksi probabilitas kejadian suatu peristiwa (respon biner). Dua pendekatan utama dalam membangun model adalah pendekatan Confirmatory dan Exploratory

1. Confirmatory (Pendekatan Konfirmator)

Pendekatan ini digunakan ketika peneliti telah memiliki teori atau hipotesis yang jelas mengenai efek atau hubungan antara variabel prediktor dan respon.

Ciri-ciri:

• Model dibangun berdasarkan teori atau hasil penelitian sebelumnya.

• Tujuan utamanya adalah menguji apakah efek tersebut benar-benar signifkan, bukan sekadar mencari model terbaik.

• Peneliti biasanya menyusun model penuh terlebih dahulu, lalu menguji apakah penambahan atau pengurangan suatu efek (misalnya, interaksi) memberikan peningkatan model secara signifkan.

• Uji signifkansi dilakukan dengan membandingkan model dengan efek tertentu dan model tanpa efek tersebut, misalnya dengan Likelihood Ratio Test.

Contoh penggunaan: Misalnya, teori menyatakan bahwa faktor x1 dan x2 mempengaruhi probabilitas seseorang membeli produk. Maka model logistik dibangun langsung dengan x1 dan x2, lalu diuji apakah kontribusi x2 benar-benar signifikan.

2. Exploratory (Pendekatan Eksploratori)

Pendekatan ini digunakan ketika peneliti belum memiliki teori yang pasti atau ingin mengeksplorasi hubungan potensial antar variabel.

Ciri-ciri:

• Model dibangun secara bertahap dengan tujuan menemukan kombinasi prediktor terbaik.

• Pemilihan variabel dilakukan berdasarkan kriteria statistik, seperti AIC, deviance, atau log-likelihood. Proses seleksi dilakukan melalui:

• Forward Selection: Mulai dari model kosong, satu per satu variabel dimasukkan jika signifkan.

• Backward Elimination: Mulai dari model penuh, variabel yang tidak signifkan dikeluarkan.

• Stepwise Selection: Gabungan dari keduanya, variabel dapat masuk dan keluar secara dinamis.

Tujuan:

Menemukan model yang parsimonious, yaitu cukup sederhana namun memiliki performa prediksi yang baik. Pemilihan antara pendekatan Confrmatory dan Exploratory bergantung pada tujuan penelitian. Jika ingin menguji hipotesis tertentu, gunakan pendekatan Confrmatory. Jika ingin menemukan model terbaik berdasarkan data, gunakan pendekatan Exploratory. Dalam praktiknya, kedua pendekatan ini sering digunakan secara komplementer: teori digunakan sebagai dasar, dan seleksi eksploratori dilakukan untuk menyempurnakan model.

Simulasi Data

Misalkan sebuah perusahaan e-commerce ingin memprediksi apakah pengunjung situs web akan membeli produk atau tidak berdasarkan karakteristik perilaku mereka saat mengakses situs. Adapun tujuan dari analisis ini adalah membangun model klasifikasi untuk memprediksi keputusan pembelian berdasarkan Waktu yang dihabiskan di halaman produk (dalam menit), Apakah pengguna berasal dari iklan media sosial (1 = ya, 0 = tidak), dan Skor ketertarikan produk yang dihitung dari klik & scroll (nilai normal).

Struktur variabelnya yaitu:

y: Faktor biner (1 = beli, 0 = tidak beli)

x1: Numerik (waktu dalam menit)

x2: Biner (asal dari iklan)

x3: Numerik (skor minat)

library(knitr)
library(dplyr)
library(ggplot2)
library(MASS)
## Warning: package 'MASS' was built under R version 4.4.3
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
library(caret)
## Warning: package 'caret' was built under R version 4.4.3
## Loading required package: lattice
library(pROC)
## Warning: package 'pROC' was built under R version 4.4.3
## Type 'citation("pROC")' for a citation.
## 
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
## 
##     cov, smooth, var
library(DescTools)
## Warning: package 'DescTools' was built under R version 4.4.3
## 
## Attaching package: 'DescTools'
## The following objects are masked from 'package:caret':
## 
##     MAE, RMSE
set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2)       # waktu di halaman
x2 <- rbinom(n, 1, 0.4)               # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1)       # skor ketertarikan

# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3

# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))

# Simulasi variabel respon
y <- rbinom(n, 1, p)

# Buat data frame
df <- data.frame(
  y = as.factor(y),
  x1 = x1,
  x2 = x2,
  x3 = x3
)

head(df)

Pemilihan Model

Model Full

model_full <- glm(y ~ x1 + x2 + x3, data = df, family = binomial)
summary(model_full)
## 
## Call:
## glm(formula = y ~ x1 + x2 + x3, family = binomial, data = df)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)   
## (Intercept)  -0.2815     0.7045  -0.400  0.68947   
## x1            0.5126     0.1597   3.210  0.00133 **
## x2            0.9536     0.5715   1.669  0.09518 . 
## x3            0.7189     0.2530   2.842  0.00448 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 134.37  on 199  degrees of freedom
## Residual deviance: 112.44  on 196  degrees of freedom
## AIC: 120.44
## 
## Number of Fisher Scoring iterations: 6

Metode Stepwise: Forward, Backward, dan Kedua Arah

null_model <- glm(y ~ 1, data = df, family = binomial)
step_forward <- step(null_model, direction = "forward", scope = formula(model_full), trace = FALSE)
step_backward <- step(model_full, direction = "backward", trace = FALSE)
step_both <- step(null_model, direction = "both", scope = formula(model_full), trace = FALSE)
AIC(model_full, step_forward, step_backward, step_both)

Evaluasi Model: ROC dan AUC

pred_prob <- predict(step_both, type = "response")
roc_obj <- roc(df$y, pred_prob)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
plot(roc_obj, main = "Kurva ROC", col = "purple")

auc(roc_obj)
## Area under the curve: 0.7768

Pseudo R Squared

PseudoR2(step_both, which = c("CoxSnell", "Nagelkerke", "McFadden"))
##   CoxSnell Nagelkerke   McFadden 
##  0.1038777  0.2123225  0.1632447

Tabel Klasifikasi dan Evaluasi

pred_class <- ifelse(pred_prob >= 0.5, 1, 0)
conf_matrix <- confusionMatrix(factor(pred_class), df$y, positive = "1")
conf_matrix
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   0   1
##          0   2   1
##          1  19 178
##                                           
##                Accuracy : 0.9             
##                  95% CI : (0.8498, 0.9378)
##     No Information Rate : 0.895           
##     P-Value [Acc > NIR] : 0.4661224       
##                                           
##                   Kappa : 0.1442          
##                                           
##  Mcnemar's Test P-Value : 0.0001439       
##                                           
##             Sensitivity : 0.99441         
##             Specificity : 0.09524         
##          Pos Pred Value : 0.90355         
##          Neg Pred Value : 0.66667         
##              Prevalence : 0.89500         
##          Detection Rate : 0.89000         
##    Detection Prevalence : 0.98500         
##       Balanced Accuracy : 0.54483         
##                                           
##        'Positive' Class : 1               
## 
conf_matrix$byClass[c("Sensitivity", "Specificity")]
## Sensitivity Specificity 
##   0.9944134   0.0952381

Metode Perbandingan Model dalam Regresi Logistik

library(MASS)
library(broom)
library(DescTools)

Simulasi Data

set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2)       # waktu di halaman
x2 <- rbinom(n, 1, 0.4)               # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1)       # skor ketertarikan

# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3

# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))

# Simulasi variabel respon
y <- rbinom(n, 1, p)

# Buat data frame
df <- data.frame(
  y = as.factor(y),
  x1 = x1,
  x2 = x2,
  x3 = x3
)

head(df)

Pembuatan Model

model1 <- glm(y ~ x1, data = df, family = binomial)
model2 <- glm(y ~ x1 + x2, data = df, family = binomial)
model3 <- glm(y ~ x1 + x2 + x3, data = df, family = binomial)

Perbandingan AIC dan Deviance

model_comp <- data.frame(
Model = c("Model 1", "Model 2", "Model 3"),
AIC = c(AIC(model1), AIC(model2), AIC(model3)),
Deviance = c(deviance(model1), deviance(model2), deviance(model3))
)
model_comp

Likelihood-Ratio Test

anova(model1, model2, test = "LRT")
anova(model2, model3, test = "LRT")

Prinsip Parsimony

Model yang kompleks sering memiliki AIC dan deviance yang lebih kecil, namun Model sederhana lebih mudah diinterpretasikan.

  • Jika penurunan AIC tidak signifikan, pilih model lebih sederhana.

  • Parsimony mencegah overfitting.

Rumus dan Penjelasan

Rumus AIC

\[ \text{AIC} = -2(\log L - k) = -2 \log L + 2k \]

Penjelasan:
AIC adalah ukuran untuk menilai model berdasarkan kombinasi antara goodness-of-fit (melalui log-likelihood) dan kompleksitas (melalui jumlah parameter \(k\)). Semakin kecil AIC, semakin baik model tersebut secara keseluruhan karena AIC menghukum model yang terlalu kompleks meskipun memiliki likelihood tinggi.

Rumus Deviance

\[ D = -2 \left[\log L(\text{model}) - \log L(\text{model saturasi})\right] \]

Penjelasan: Deviance mengukur seberapa jauh model saat ini dibandingkan dengan model sempurna (model saturasi). Nilai deviance yang kecil menu

Rumus Likelihood-Ratio

\[ G^2 = -2 (\log L_0 - \log L_1) \]

Penjelasan:
Statistik Likelihood Ratio digunakan untuk menguji apakah penambahan variabel dalam model secara signifikan meningkatkan kecocokan model. Jika \(G^2\) besar dan p-value kecil, maka model kompleks lebih baik dari model sederhana secara statistik.

Evaluasi Tabel Klasifkasi dan Akurasi Model

pred_prob <- predict(model3, type = "response")
pred_class <- factor(ifelse(pred_prob >= 0.5, 1, 0))
conf_matrix <- confusionMatrix(pred_class, df$y, positive = "1")
conf_matrix
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   0   1
##          0   2   1
##          1  19 178
##                                           
##                Accuracy : 0.9             
##                  95% CI : (0.8498, 0.9378)
##     No Information Rate : 0.895           
##     P-Value [Acc > NIR] : 0.4661224       
##                                           
##                   Kappa : 0.1442          
##                                           
##  Mcnemar's Test P-Value : 0.0001439       
##                                           
##             Sensitivity : 0.99441         
##             Specificity : 0.09524         
##          Pos Pred Value : 0.90355         
##          Neg Pred Value : 0.66667         
##              Prevalence : 0.89500         
##          Detection Rate : 0.89000         
##    Detection Prevalence : 0.98500         
##       Balanced Accuracy : 0.54483         
##                                           
##        'Positive' Class : 1               
## 

Sensitivitas dan Spesifisitas

  • Sensitivitas: Kemampuan model mendeteksi kelas positif secara benar (True Positive Rate)

    \[ \text{Sensitivity} = \frac{\text{TP}}{\text{TP} + \text{FN}} \]

  • Spesifisitas: Kemampuan model mendeteksi kelas negatif secara benar (True Negative Rate)

    \[ \text{Specificity} = \frac{\text{TN}}{\text{TN} + \text{FP}} \]

    conf_matrix$byClass[c("Sensitivity", "Specificity")]
    ## Sensitivity Specificity 
    ##   0.9944134   0.0952381

Kesimpulan

• Deviance yang kecil menunjukkan kecocokan model yang lebih baik.

• AIC yang rendah menunjukkan keseimbangan antara kecocokan dan kompleksitas.

• Likelihood Ratio Test mengevaluasi apakah model kompleks secara signifkan lebih baik.

• Tabel klasifkasi membantu menilai kinerja prediksi aktual vs prediksi model.

• Prinsip Parsimony mengutamakan model sederhana jika performanya mirip.

Detail ROC & Penjelasan Kurva ROC (Receiver Operating Characteristic)

Kurva ROC adalah alat visual yang digunakan untuk mengevaluasi performa model klasifkasi biner. Kurva ini menunjukkan trade-of antara True Positive Rate (Sensitivity) dan False Positive Rate (1 - Specificity) pada berbagai threshold klasifkasi.

  1. Definisi

• Sumbu Y: Sensitivity = True Positive Rate = TP / (TP + FN)

• Sumbu X: 1 - Specifcity = False Positive Rate = FP / (FP + TN)

• Garis diagonal (dari kiri bawah ke kanan atas) menunjukkan performa acak (random guess).

• Kurva yang mendekati pojok kiri atas menunjukkan performa klasifkasi yang lebih baik.

  1. Cut-off dan Pergerakan Kurva

Saat cut-off menurun, model mengklasifkasikan lebih banyak pengamatan sebagai positif:

• Sensitivitas naik

• Spesifisitas turun

Saat cut-off naik, model menjadi lebih konservatif:

• Sensitivitas turun

• Spesifisitas naik

  1. Kurva ROC Ideal

Kurva ideal memiliki bentuk:

• Naik tajam secara vertikal hingga mencapai sensitivitas = 1

• Lalu bergerak secara horizontal menuju 1 - specifcity = 1

• Area under the curve (AUC) mendekati 1

  1. Interpretasi Luas Area (AUC)

• AUC = 0.5: model tidak lebih baik dari tebak acak

• AUC > 0.7: model cukup baik

• AUC > 0.9: model sangat baik

• AUC dikenal juga sebagai concordance index, yaitu probabilitas bahwa model memberikan nilai skor probabilitas yang lebih tinggi untuk kasus positif daripada kasus negatif.

  1. Kegunaan Kurva ROC

• Untuk membandingkan performa beberapa model klasifkasi

• Untuk memilih threshold (cut-of) optimal berdasarkan kebutuhan aplikasi (misalnya: lebih penting menghindari false negative atau false positive)

  1. Visualisasi dalam R Kurva ROC dapat dibuat menggunakan package pROC:
library(pROC)
set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2)       # waktu di halaman
x2 <- rbinom(n, 1, 0.4)               # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1)       # skor ketertarikan

# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3

# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))

# Simulasi variabel respon
y <- rbinom(n, 1, p)

data <- data.frame(y = as.factor(y), x1, x2, x3)
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
pred <- predict(model, type = "response")
roc_obj <- roc(data$y, pred)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
plot(roc_obj)

auc(roc_obj)
## Area under the curve: 0.7768
  1. Simulasi Pemilihan Threshold Optimal

Untuk memilih threshold terbaik, kita bisa mengevaluasi sensitivitas dan spesifsitas pada berbagai cut-off.

thresholds <- seq(0.1, 0.9, by = 0.05)
results <- data.frame(Threshold = thresholds)
true_y <- factor(data$y, levels = c(0, 1))
results$Sensitivity <- sapply(thresholds, function(t) {
  pred_class <- factor(ifelse(pred >= t, 1, 0), levels = c(0, 1))
  cm <- table(Pred = pred_class, Obs = true_y)

  TP <- ifelse("1" %in% rownames(cm) && "1" %in% colnames(cm), cm["1", "1"], 0)
  FN <- ifelse("0" %in% rownames(cm) && "1" %in% colnames(cm), cm["0", "1"], 0)

  if ((TP + FN) == 0) return(NA)
  TP / (TP + FN)
})

results$Specificity <- sapply(thresholds, function(t) {
  pred_class <- factor(ifelse(pred >= t, 1, 0), levels = c(0, 1))
  cm <- table(Pred = pred_class, Obs = true_y)

  TN <- ifelse("0" %in% rownames(cm) && "0" %in% colnames(cm), cm["0", "0"], 0)
  FP <- ifelse("1" %in% rownames(cm) && "0" %in% colnames(cm), cm["1", "0"], 0)

  if ((TN + FP) == 0) return(NA)
  TN / (TN + FP)
})

print(results)
##    Threshold Sensitivity Specificity
## 1       0.10   1.0000000  0.00000000
## 2       0.15   1.0000000  0.00000000
## 3       0.20   1.0000000  0.00000000
## 4       0.25   1.0000000  0.00000000
## 5       0.30   1.0000000  0.00000000
## 6       0.35   1.0000000  0.04761905
## 7       0.40   1.0000000  0.04761905
## 8       0.45   0.9944134  0.09523810
## 9       0.50   0.9944134  0.09523810
## 10      0.55   0.9888268  0.14285714
## 11      0.60   0.9832402  0.19047619
## 12      0.65   0.9832402  0.28571429
## 13      0.70   0.9664804  0.33333333
## 14      0.75   0.9553073  0.33333333
## 15      0.80   0.8938547  0.38095238
## 16      0.85   0.8212291  0.61904762
## 17      0.90   0.6480447  0.71428571

Cut-off optimal bisa dipilih berdasarkan:

Maksimum dari Sensitivity + Specificity

Atau mempertimbangkan trade-of sesuai tujuan aplikasi (misalnya: jika False Negative harus dihindari, maka prioritaskan sensitivitas tinggi)

  1. Catatan

ROC cocok saat proporsi kelas seimbang

Untuk data dengan kelas tidak seimbang, precision-recall curve bisa lebih informatif

Precision-Recall Curve (PR Curve)

Penjelasan Precision-Recall Curve

Kurva Precision-Recall (PR) adalah alat evaluasi performa model klasifikasi, khususnya sangat berguna saat bekerja dengan data yang tidak seimbang (class imbalance).

1. Definisi

  • Precision (Presisi): Proporsi prediksi positif yang benar-benar positif

    \[ \text{Precision} = \frac{\text{TP}}{\text{TP} + \text{FP}} \]

  • Recall (Sensitivitas): Proporsi kasus positif yang berhasil diprediksi positif

    \[ \text{Recall} = \frac{\text{TP}}{\text{TP} + \text{FN}} \]

2. Interpretasi

  • PR Curve menunjukkan bagaimana presisi berubah saat recall meningkat.
  • Idealnya, kita ingin nilai presisi dan recall keduanya tinggi, tetapi biasanya ada trade-off.
  • Model dengan performa baik memiliki PR Curve yang melengkung ke pojok kanan atas.

3. Area Under PR Curve

  • Luas kurva (AUPRC) mendekati 1 berarti model sangat baik.
  • Baseline AUPRC = prevalensi kelas positif dalam data.

4. PR Curve vs ROC Curve

Aspek ROC Curve Precision-Recall Curve
Fokus Semua kelas Kelas positif saja
Kuat di Data seimbang Data tidak seimbang
Sumbu Y Sensitivitas (Recall) Precision
Sumbu X 1 - Spesifisitas Recall

5. Visualisasi PR Curve di R

library(PRROC)
## Warning: package 'PRROC' was built under R version 4.4.3
## Loading required package: rlang
set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2)       # waktu di halaman
x2 <- rbinom(n, 1, 0.4)               # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1)       # skor ketertarikan

# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3

# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))

# Simulasi variabel respon
y <- rbinom(n, 1, p)
data <- data.frame(y = y, x1, x2, x3)
head(data)
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
prob <- predict(model, type = "response")
pr <- pr.curve(scores.class0 = prob[data$y == 1],
scores.class1 = prob[data$y == 0],
curve = TRUE)
plot(pr)

6. Catatan

• PR Curve sangat informatif untuk aplikasi seperti deteksi penipuan atau diagnosis penyakit langka.

• PR Curve digunakan saat:

– Kelas positif jauh lebih jarang daripada kelas negatif

– Tujuan aplikasi lebih mementingkan presisi terhadap kelas minoritas

Pseudo R-squared pada Regresi Logistik

Simulasi Data

set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2)       # waktu di halaman
x2 <- rbinom(n, 1, 0.4)               # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1)       # skor ketertarikan

# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3

# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))

# Simulasi variabel respon
y <- rbinom(n, 1, p)
data <- data.frame(y = as.factor(y), x1, x2, x3)
head(data)

Model Logistik dan Null Model

model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
model_null <- glm(y ~ 1, data = data, family = binomial)

Likelihood dan Rumus \[ R^2_{\text{Cox and Snell}} = 1 - \left( \frac{L_0}{L_M} \right)^{2/n} \]

\[ R^2_{\text{McFadden}} = 1 - \frac{\log L_M}{\log L_0} \]

Dengan:

  • \(L_0\): likelihood model null (tanpa prediktor)
  • \(L_M\): likelihood model penuh

Perhitungan Manual R-squared

logL0 <- logLik(model_null)
logLM <- logLik(model)
L0 <- exp(logL0)
LM <- exp(logLM)
n <- nobs(model)
cox_snell <- 1 - (L0 / LM)^(2 / n)
mcfadden <- 1 - (as.numeric(logLM) / as.numeric(logL0))
r2 <- data.frame(
R2_Cox_Snell = cox_snell,
R2_McFadden = mcfadden
)
r2

Perhitungan Otomatis dengan Package Tambahan

Menggunakan DescTools

if (!require(DescTools)) install.packages("DescTools"); library(DescTools)
PseudoR2(model, which = "all")
##        McFadden     McFaddenAdj        CoxSnell      Nagelkerke   AldrichNelson 
##      0.16324473      0.10370891      0.10387765      0.21232254      0.09883794 
## VeallZimmermann           Efron McKelveyZavoina            Tjur             AIC 
##      0.24594790      0.14859709      0.31510399      0.14005131    120.43722088 
##             BIC          logLik         logLik0              G2 
##    133.63049035    -56.21861044    -67.18644287     21.93566486

Interpretasi

Nilai McFadden R² sebesar 0.163 menunjukkan bahwa model memiliki kecocokan yang moderat.

Nilai Cox & Snell R² sebesar 0.1039 bersifat lebih konservatif dan memang tidak bisa mencapai 1 penuh karena keterbatasan metode ini. Nilai ini menunjukkan bahwa sekitar 10.4% variasi dalam data dijelaskan oleh model.

Nilai Nagelkerke R² sebesar 0.2123 merupakan versi koreksi dari Cox & Snell agar dapat mencapai maksimum 1. Dengan nilai > 0.2, ini memberikan indikasi bahwa model memiliki kecocokan yang layak

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 satu dari tiga jenis buah favorit: Apel (A), Jeruk (J), dan Pisang (P).

Hasil survei:

  • Apel: 4 orang
  • Jeruk: 3 orang
  • Pisang: 3 orang

Probabilitas teoretik preferensi:

  • \(p_A = 0.3\)
  • \(p_J = 0.4\)
  • \(p_P = 0.3\)

Pertanyaannya: Berapa peluang bahwa dalam 10 orang akan ada 4 yang memilih Apel, 3 memilih Jeruk, dan 3 memilih Pisang?

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 = 10\), \(x_1 = 4\), \(x_2 = 3\), \(x_3 = 3\)
  • \(p_1 = 0.3\), \(p_2 = 0.4\), \(p_3 = 0.3\)

Perhitungan Manual di R

n <- 10
x <- c(4, 3, 3)
p <- c(0.3, 0.4, 0.3)
# Hitung komponen-koefisien
faktorial_total <- factorial(n)
faktorial_x <- prod(factorial(x))
koefisien <- faktorial_total / faktorial_x
# Hitung peluang
peluang <- koefisien * prod(p^x)
peluang
## [1] 0.05878656

Probabilitas bahwa 4 orang memilih Apel, 3 Jeruk, dan 3 Pisang (dengan proporsi preferensi 0.3, 0.4, dan 0.3) adalah 0.058787. Distribusi multinomial digunakan untuk menghitung peluang dalam percobaan dengan beberapa kategori hasil. Rumus dasarnya merupakan generalisasi dari binomial untuk lebih dari dua kategori.

Multinomial Logistic 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_{0} + \beta_{1j}x_{1} + \cdots + \beta_{pj}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( \frac{\pi_1 / \pi_3}{\pi_2 / \pi_3} \right) = \log \left( \frac{\pi_1}{\pi_3} \right) - \log \left( \frac{\pi_2}{\pi_3} \right) \] \[ = (\alpha_1 + \beta_1 x) - (\alpha_2 + \beta_2 x) = (\alpha_1 - \alpha_2) + (\beta_1 - \beta_2)x \] 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 ditentukan dengan relevel().

Estimasi Parameter

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

Log-likelihood:

\[ l(\beta) = \sum_{i=1}^n \sum_{j=1}^K y_{ij} \log (\pi_{ij}) \]

dengan \(\pi_{ij} = P(Y = j | x_i)\) dan \(y_{ij} = 1\) jika \(Y_i = j\), 0 jika tidak.

Contoh Kasus

Sebuah institusi pendidikan ingin memahami faktor-faktor yang memengaruhi preferensi gaya belajar siswa, yaitu Visual, Auditory, atau Kinesthetic. Sebanyak 150 siswa disurvei, dan dikumpulkan data sebagai berikut:

  • LearningStyle: Gaya belajar yang dipilih (Visual, Auditory, Kinesthetic)

  • GradeLevel: Tingkatan kelas siswa (Junior, Senior, Freshman)

  • StudyHours: Rata-rata jam belajar mandiri per minggu

  • Extracurriculars: Jumlah kegiatan ekstrakurikuler yang diikuti

Tujuan analisis ini adalah untuk mengetahui bagaimana GradeLevel, StudyHours, dan Extracurriculars memengaruhi preferensi gaya belajar siswa.

Simulasi Data

set.seed(123)
n <- 150

GradeLevel <- sample(c("Freshman", "Junior", "Senior"), n, replace = TRUE)
StudyHours <- round(rnorm(n, mean = 10, sd = 3))
Extracurriculars <- sample(0:5, n, replace = TRUE)

# Simulasikan gaya belajar berdasarkan grade level
LearningStyle <- sapply(GradeLevel, function(level) {
  if (level == "Freshman") {
    sample(c("Visual", "Auditory", "Kinesthetic"), size = 1, prob = c(0.5, 0.3, 0.2))
  } else if (level == "Junior") {
    sample(c("Visual", "Auditory", "Kinesthetic"), size = 1, prob = c(0.3, 0.4, 0.3))
  } else {
    sample(c("Visual", "Auditory", "Kinesthetic"), size = 1, prob = c(0.2, 0.3, 0.5))
  }
})

# Buat data frame
df <- data.frame(
  LearningStyle = factor(LearningStyle),
  GradeLevel = factor(GradeLevel),
  StudyHours,
  Extracurriculars
)

# Tetapkan baseline untuk analisis (misal: Visual sebagai baseline)
df$LearningStyle <- relevel(df$LearningStyle, ref = "Visual")

head(df)

Estimasi Model

library(nnet)
model_mnlogit <- multinom(LearningStyle ~ GradeLevel + StudyHours + Extracurriculars, data = df)
## # weights:  18 (10 variable)
## initial  value 164.791843 
## iter  10 value 150.818303
## final  value 150.779592 
## converged
summary(model_mnlogit)
## Call:
## multinom(formula = LearningStyle ~ GradeLevel + StudyHours + 
##     Extracurriculars, data = df)
## 
## Coefficients:
##             (Intercept) GradeLevelJunior GradeLevelSenior StudyHours
## Auditory      -1.313043         1.176353        0.9339879 0.09852580
## Kinesthetic   -1.977086         1.818910        2.1623436 0.03622559
##             Extracurriculars
## Auditory           0.0101990
## Kinesthetic        0.2123591
## 
## Std. Errors:
##             (Intercept) GradeLevelJunior GradeLevelSenior StudyHours
## Auditory      0.8645914        0.5164149        0.5376260 0.07033911
## Kinesthetic   0.9351475        0.5908354        0.5904723 0.07194458
##             Extracurriculars
## Auditory           0.1330148
## Kinesthetic        0.1377687
## 
## Residual Deviance: 301.5592 
## AIC: 321.5592

Nilai P-Value dan Interpretasi

z <- summary(model_mnlogit)$coefficients / summary(model_mnlogit)$standard.errors
pval <- 2 * (1 - pnorm(abs(z)))
round(pval, 4)
##             (Intercept) GradeLevelJunior GradeLevelSenior StudyHours
## Auditory         0.1288           0.0227           0.0823     0.1613
## Kinesthetic      0.0345           0.0021           0.0003     0.6146
##             Extracurriculars
## Auditory              0.9389
## Kinesthetic           0.1232

Interpretasi:

• Koefisien untuk kategori “Auditory” dan “Kinesthetic” dibandingkan dengan baseline “Visual”

• Nilai p-value kecil (<0.05) menunjukkan variabel tersebut signifkan memengaruhi preferensi gaya belajar.

Prediksi dan Validasi

df$Predicted <- predict(model_mnlogit)
table(Predicted = df$Predicted, Actual = df$LearningStyle)
##              Actual
## Predicted     Visual Auditory Kinesthetic
##   Visual          18       12           4
##   Auditory        10       26          17
##   Kinesthetic     11       18          34

Kesimpulan

Model regresi logistik multinomial berhasil digunakan untuk:

• Menganalisis hubungan antara atribut sisawa dan preferensi gaya belajar

• Mengetahui faktor signifkan yang memengaruhi pilihan

• Memungkinkan prediksi gaya belajar yang dipilih oleh siswa berdasarkan karakteristiknya

Contoh Kasus 2 di R

Kita gunakan dataset iris untuk memodelkan jenis spesies bunga berdasarkan lebar dan panjang kelopak.

data(iris)
iris <- iris %>% mutate(Species = relevel(Species, ref = "setosa"))
model <- multinom(Species ~ Petal.Length + Petal.Width, data = iris)
## # weights:  12 (6 variable)
## initial  value 164.791843 
## iter  10 value 12.657828
## iter  20 value 10.374056
## iter  30 value 10.330881
## iter  40 value 10.306926
## iter  50 value 10.300057
## iter  60 value 10.296452
## iter  70 value 10.294046
## iter  80 value 10.292029
## iter  90 value 10.291154
## iter 100 value 10.289505
## final  value 10.289505 
## stopped after 100 iterations
summary(model)
## Call:
## multinom(formula = Species ~ Petal.Length + Petal.Width, data = iris)
## 
## Coefficients:
##            (Intercept) Petal.Length Petal.Width
## versicolor   -22.79944      6.92122    7.878496
## virginica    -67.82521     12.64721   18.261016
## 
## Std. Errors:
##            (Intercept) Petal.Length Petal.Width
## versicolor     44.3859     37.58715    81.00888
## virginica      46.3939     37.65702    81.09482
## 
## Residual Deviance: 20.57901 
## AIC: 32.57901

Interpretasi Koefisien

z <- summary(model)$coefficients / summary(model)$standard.errors
p_values <- 2 * (1 - pnorm(abs(z)))
round(p_values, 4)
##            (Intercept) Petal.Length Petal.Width
## versicolor      0.6075       0.8539      0.9225
## virginica       0.1438       0.7370      0.8218

Nilai p-value menunjukkan apakah variabel prediktor berpengaruh signifkan terhadap log odds dibanding baseline category.

Prediksi dan Visualisasi

iris$predicted <- predict(model, newdata = iris)
table(Predicted = iris$predicted, Actual = iris$Species)
##             Actual
## Predicted    setosa versicolor virginica
##   setosa         50          0         0
##   versicolor      0         47         3
##   virginica       0          3        47
ggplot(iris, aes(x = Petal.Length, y = Petal.Width, color = predicted)) +
geom_point(size = 2) +
labs(title = "Multinomial Logistic Regression Predictions",
x = "Petal Length", y = "Petal Width") +
theme_minimal()

Kesimpulan Multinomial logistic regression efektif digunakan untuk klasifikasi kategori nominal. Model ini memberikan estimasi log-odds terhadap baseline dan prediksi kategori baru berdasarkan kombinasi prediktor.

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 ditulis:

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

Contoh Data: Kepuasan Pelanggan

Misalkan seorang dosen ingin mengetahui apakah semakin tinggi intensitas belajar (jam belajar per minggu) mahasiswa berhubungan dengan tingkat pemahaman mereka terhadap materi.

Variabel yang digunakan:

  • understanding (ordinal):

1: Rendah

2: Sedang

3: Tinggi

  • study_hours (numerik): jumlah jam belajar per minggu
set.seed(123)
n <- 200
study_hours <- round(runif(n, 1, 20))  # Jam belajar antara 1 hingga 20 jam
understanding <- cut(3 + 0.3*study_hours + rnorm(n),
                     breaks = c(-Inf, 6.5, 9.5, Inf),
                     labels = c("Rendah", "Sedang", "Tinggi"),
                     ordered_result = TRUE)
df <- data.frame(understanding, study_hours)
head(df)

Estimasi Model Ordinal

model_ord <- polr(understanding ~ study_hours, data = df, Hess = TRUE)
summary(model_ord)
## Call:
## polr(formula = understanding ~ study_hours, data = df, Hess = TRUE)
## 
## Coefficients:
##              Value Std. Error t value
## study_hours 0.5241    0.06772    7.74
## 
## Intercepts:
##               Value   Std. Error t value
## Rendah|Sedang  6.1901  0.8153     7.5923
## Sedang|Tinggi 11.9717  1.3113     9.1299
## 
## Residual Deviance: 162.6821 
## AIC: 168.6821

Nilai P-Value

(ctable <- coef(summary(model_ord)))
##                    Value Std. Error  t value
## study_hours    0.5240819 0.06771506 7.739519
## Rendah|Sedang  6.1900621 0.81531258 7.592256
## Sedang|Tinggi 11.9717293 1.31126286 9.129923
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
## study_hours    0.5240819 0.06771506 7.739519       0
## Rendah|Sedang  6.1900621 0.81531258 7.592256       0
## Sedang|Tinggi 11.9717293 1.31126286 9.129923       0

Prediksi Probabilitas

newdata <- data.frame(study_hours = 5:9)
predict(model_ord, newdata = newdata, type = "probs")
##      Rendah     Sedang       Tinggi
## 1 0.9726059 0.02730723 0.0000868432
## 2 0.9545909 0.04526244 0.0001466614
## 3 0.9256346 0.07411775 0.0002476725
## 4 0.8805245 0.11905726 0.0004182246
## 5 0.8135622 0.18573168 0.0007061392

Goodness-of-Fit dan Proportional Odds

Model cumulative logit mengasumsikan efek prediktor sama untuk setiap cutoff. 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 devian 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).

Asumsi paralelisme menyatakan bahwa koefisien regresi (\(\beta\)) 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, 2, \dots, c - 1\):

  • Hanya intercept (\(\alpha_j\)) yang berbeda-beda.
  • Koefisien \(\beta\) tetap sama untuk semua logit fungsi 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 digunakan 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 (pada regresi ordinal)

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 kategorikal memegang peranan penting dalam statistika terapan karena banyak kejadian di dunia nyata menghasilkan data dalam kategori, seperti jenis kelamin, status pekerjaan, tingkat pendidikan, preferensi konsumen, atau diagnosis medis. Data yang berbentuk kategori ini biasanya dianalisis dengan menggunakan tabel kontingensi, model log-linier, dan model regresi logistik. Setiap metode memiliki keuntungan dan kerugian yang bervariasi tergantung pada tujuan analisis serta struktur data.

Tabel kontingensi berfungsi sebagai langkah awal untuk mengeksplorasi hubungan antara dua atau lebih variabel kategorikal. Contohnya, dalam penelitian mengenai dampak obat terhadap serangan jantung, tabel kontingensi dapat menunjukkan jumlah pasien yang mengalami atau tidak mengalami serangan jantung berdasarkan obat yang digunakan. Tabel ini berperan dalam mengidentifikasi pola-pola awal serta menghitung ukuran asosiasi seperti odds ratio, risk ratio, dan statistik chi-square untuk menguji independensi antara variabel.

Namun, untuk menciptakan model statistik yang dapat mengontrol efek dari berbagai variabel dan interaksinya secara bersamaan, model log-linier menjadi sangat bermanfaat. Model log-linier merupakan tipe khusus dari Generalized Linear Model (GLM) yang diterapkan pada frekuensi sel dalam tabel kontingensi dan mengasumsikan distribusi Poisson. Berbeda dengan regresi logistik, model log-linier tidak menentukan variabel mana yang dependen dan mana yang independen, melainkan memperlakukannya secara simetris. Model ini lebih sesuai bila tujuan analisis adalah untuk memahami struktur asosiasi atau independensi antar variabel, alih-alih untuk tujuan prediksi.

Struktur model log-linier ditentukan oleh efek utama dari masing-masing variabel serta interaksi di antara mereka. Sebagai contoh, dalam tabel kontingensi tiga arah (seperti: jenis kelamin, status merokok, dan penyakit paru-paru), model ini mampu mengidentifikasi apakah interaksi antara dua variabel cukup untuk menjelaskan data, atau apakah interaksi tiga arah diperlukan untuk memahami struktur asosiasi yang ada. Penyesuaian model bisa dilakukan melalui metode uji rasio kemungkinan untuk membandingkan model yang lebih sederhana dengan yang lebih rumit.

Sebaliknya, regresi logistik adalah metode yang paling umum dipakai ketika satu variabel kategorikal secara spesifik dianggap sebagai variabel dependen (misalnya, kejadian penyakit: ya/tidak) dan satu atau lebih variabel kategorik maupun numerik berfungsi sebagai prediktor. Model ini mengkaji logit dari probabilitas kejadian (atau log odds), dan sangat berguna dalam penelitian observasional serta eksperimental untuk menjelaskan atau meramalkan kemungkinan suatu hasil. Regresi logistik juga memiliki variasi untuk outcome kategorik yang lebih dari dua kelas, termasuk regresi logistik multinomial dan regresi logistik ordinal.

Oleh karena itu, meskipun ketiga metode ini beroperasi pada data kategorikal, tabel kontingensi bersifat deskriptif, model log-linier bersifat eksploratif dalam hubungan simetris, sedangkan regresi logistik adalah alat prediktif untuk outcome kategorikal. Pemilihan metode yang tepat sangat bergantung pada apakah analisis berfokus pada deskripsi, eksplorasi struktur, atau prediksi hasil berdasarkan variabel penjelas. Paduan ketiga metode ini sering dijumpai dalam praktik untuk memperoleh pemahaman yang lebih mendalam tentang data kategorikal yang dianalisis.

Ringkasan Dalam analisis data kategorikal, terdapat beberapa metode statistik yang umum dipakai, di antaranya: 1. Tabel Kontingensi: menampilkan frekuensi gabungan dari dua atau lebih variabel kategorikal. 2. Model Log-linier: digunakan untuk memodelkan struktur asosiasi di dalam tabel kontingensi tanpa menganggap adanya variabel dependen. 3. Model Regresi Logistik: digunakan untuk memodelkan probabilitas dari kategori variabel dependen berdasarkan variabel independen.

Meskipun ketiga metode ini dapat diaplikasikan pada data kategorikal, pendekatan dan cara interpretasinya sangatlah berbeda.

Tabel Kontingensi

Tabel kontingensi menyajikan jumlah frekuensi dari kombinasi kategori antar variabel.

Contoh tabel 2x2:

# Membuat matriks dari data 
table_data <- matrix(c(28, 656, 18, 658), 
                     nrow = 2, 
                     byrow = TRUE,
                     dimnames = list(Obat = c("Placebo", "Aspirin"),
                                     Infeksi = c("Ya", "Tidak")))
table_data
##          Infeksi
## Obat      Ya Tidak
##   Placebo 28   656
##   Aspirin 18   658

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}) = \mu + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]

  • 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 * Infeksi, data = table_data)
## Call:
## loglm(formula = ~Obat * Infeksi, 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(
  Infeksi = c(1, 0, 1, 0),
  Obat = factor(c("Placebo", "Placebo", "Aspirin", "Aspirin")),
  Frek = c(28, 656, 18, 658)
)
head(data_glm)
model_logit <- glm(Infeksi ~ Obat, weights = Frek, family = binomial, data = data_glm)
summary(model_logit)
## 
## Call:
## glm(formula = Infeksi ~ Obat, family = binomial, data = data_glm, 
##     weights = Frek)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -3.5988     0.2389 -15.065   <2e-16 ***
## ObatPlacebo   0.4449     0.3071   1.449    0.147    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 401.99  on 3  degrees of freedom
## Residual deviance: 399.85  on 2  degrees of freedom
## AIC: 403.85
## 
## Number of Fisher Scoring iterations: 6

Perbandingan Ketiga Pendekatan

Aspek Tabel Kontingensi Model Loglinear Regresi Logistik
Tujuan Deskripsi frekuensi Deteksi asosiasi Prediksi probabilitas
Variabel dependen Tidak ada Tidak ada (intersit) Ada (eksplisit)
Distribusi Tidak diasumsikan Poisson (frekuensi sel) Binomial (probabilitas)
Bentuk Model Tidak ada GLM: log(μ) ~ efek GLM: logit(p) ~ prediktor
Cocok untuk Eksplorasi awal Tabel > 2 variabel Studi prediktif

Tabel Kontingensi dan Model Loglinier

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

# Contoh tabel 2x2
table_data <- matrix(c(28, 656, 18, 658), 
                     nrow = 2, 
                     byrow = TRUE,
                     dimnames = list(Obat = c("Placebo", "Aspirin"),
                                     Infeksi = c("Ya", "Tidak")))
table_data
##          Infeksi
## Obat      Ya Tidak
##   Placebo 28   656
##   Aspirin 18   658

Model log-linier untuk tabel I x J dapat dituliskan: \[ \log(\mu_{ij}) = \mu + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]

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(28, 656, 18, 658), 
                     nrow = 2, 
                     byrow = TRUE,
                     dimnames = list(Obat = c("Placebo", "Aspirin"),
                                     Infeksi  = c("Ya", "Tidak")))
ftable(data)
##         Infeksi  Ya Tidak
## Obat                     
## Placebo          28   656
## Aspirin          18   658

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

model_saturated <- loglm(~ Obat * Infeksi, data = data)
summary(model_saturated)
## Formula:
## ~Obat * Infeksi
## attr(,"variables")
## list(Obat, Infeksi)
## attr(,"factors")
##         Obat Infeksi Obat:Infeksi
## Obat       1       0            1
## Infeksi    0       1            1
## attr(,"term.labels")
## [1] "Obat"         "Infeksi"      "Obat:Infeksi"
## 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}) = \mu + \lambda^T_i + \lambda^B_j \] Model ini menguji hipotesis bahwa variabel X dan Y saling independen.

model_indep <- loglm(~ Obat + Infeksi, data = data)
summary(model_indep)
## Formula:
## ~Obat + Infeksi
## attr(,"variables")
## list(Obat, Infeksi)
## attr(,"factors")
##         Obat Infeksi
## Obat       1       0
## Infeksi    0       1
## attr(,"term.labels")
## [1] "Obat"    "Infeksi"
## 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.147353  1 0.1428159
## Pearson          2.129972  1 0.1444434

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.4448769

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 + Infeksi 
## Model 2:
##  ~Obat * Infeksi 
## 
##           Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1   2.147353  1                                    
## Model 2   0.000000  0   2.147353         1        0.14282
## Saturated 0.000000  0   0.000000         0        1.00000

Studi Kasus

data_survey <- matrix(
  c(25, 25, 12, 0, 1, 3),
  nrow = 2,
  byrow = TRUE,
  dimnames = list(
    Perilaku = c("Dikontrol", "Infeksi"),
    Rokok = c("0", "1-24", ">25")
  )
)
ftable(data_survey)
##           Rokok  0 1-24 >25
## Perilaku                   
## Dikontrol       25   25  12
## Infeksi          0    1   3
loglm(~ Perilaku + Rokok, data = data_survey)
## Call:
## loglm(formula = ~Perilaku + Rokok, data = data_survey)
## 
## Statistics:
##                       X^2 df   P(> X^2)
## Likelihood Ratio 6.690106  2 0.03525835
## Pearson          6.956203  2 0.03086595

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 (\(\mu_{ij}\)) dapat dinyatakan sebagai penjumlahan efek variabel dan (bila perlu) interaksinya. Untuk tabel 2x2

\[ \log(\mu_{ij}) = \lambda + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]

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{aligned} &\lambda^A_1 + \lambda^A_2 = 0 \\ &\lambda^B_1 + \lambda^B_2 = 0 \\ &\lambda^{AB}_{11} + \lambda^{AB}_{12} + \lambda^{AB}_{21} + \lambda^{AB}_{22} = 0 \end{aligned} \]

Rumus Estimasi Parameter dengan Sum-to-Zero Constraint

\[ \lambda^A_1 = \frac{1}{2} \left[ (\log \mu_{11} + \log \mu_{12}) - (\log \mu_{21} + \log \mu_{22}) \right] \]

\[ \lambda^B_1 = \frac{1}{2} \left[ (\log \mu_{11} + \log \mu_{21}) - (\log \mu_{12} + \log \mu_{22}) \right] \]

\[ \lambda^{AB}_{12} = \frac{1}{4} \left[ \log \mu_{12} - \log \mu_{11} - \log \mu_{22} + \log \mu_{21} \right] \]

Analisis Data Tabel Kontingensi 2x2

Diberikan data:

Infeksi Tidak Infeksi
Placebo 28 656
Aspirin 18 658

Bentuk Model Log-Linear

Model log-linear pada tabel 2x2:

\[ \log(\mu_{ij}) = \lambda + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]

dengan constraint sum-to-zero:

\[ \sum_i \lambda^A_i = 0,\quad \sum_j \lambda^B_j = 0,\quad \sum_{i,j} \lambda^{AB}_{ij} = 0 \]

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

Misalkan:

  • A1 = Aspirin, A2 = Placebo
  • B1 = Infeksi, B2 = Tidak Infeksi

Observasi:

  • \(n_{11} = 28\), \(n_{12} = 656\)
  • \(n_{21} = 18\), \(n_{22} = 658\)

\[ \begin{aligned} \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{aligned} \]

Dengan constraint sum-to-zero:

\[ \begin{aligned} \lambda^A_1 + \lambda^A_2 &= 0, \\ \lambda^B_1 + \lambda^B_2 &= 0, \\ \lambda^{AB}_{11} + \lambda^{AB}_{12} + \lambda^{AB}_{21} + \lambda^{AB}_{22} &= 0 \end{aligned} \]

Langkah-langkah:

1. Hitung rata-rata log frekuensi sel:

\[ \lambda = \frac{1}{4} \sum_{i=1}^{2} \sum_{j=1}^{2} \log(n_{ij}) \]

\[ = \frac{1}{4} (\log(28) + \log(656) + \log(18) + \log(658)) \]

\[ = \frac{1}{4} (3.3322 + 6.4850 + 2.8904 + 6.4887) \]

\[ = 4.7991 \]

2. Efek utama A (Jenis Obat):

\[ \lambda^A_1 = \frac{1}{2} [(\log(28) + \log(656)) - (\log(18) + \log(658))] \]

\[ = \frac{1}{2} [(3.3322 + 6.4850) - (2.8904 + 6.4887)] \]

\[ = \frac{1}{2} (0.4381) = 0.2190 \]

\[ \lambda^A_2 = -0.2190 \]

3. Efek utama B (Kategori Infeksi):

\[ \lambda^B_1 = \frac{1}{2} [(\log(28) + \log(18)) - (\log(656) + \log(658))] \]

\[ = \frac{1}{2} [(3.3322 + 2.8904) - (6.4850 + 6.4887)] \]

\[ = \frac{1}{2} (6.2226 - 12.9737) = \frac{1}{2} (-6.7511) = -3.3756 \]

\[ \lambda^B_2 = +3.3756 \]

4. Efek interaksi:

\[ \lambda^{AB}_{11} = \frac{1}{4} [\log(28) - \log(656) - \log(18) + \log(658)] \]

\[ = \frac{1}{4} [3.3322 - 6.4850 - 2.8904 + 6.4887] \]

\[ = \frac{1}{4} (0.4455) = 0.1114 \]

\[ \lambda^{AB}_{12} = -0.1114 \quad \lambda^{AB}_{21} = -0.1114 \quad \lambda^{AB}_{22} = +0.1114 \]

Ringkasan parameter:

\[ \lambda = 4.7991 \]

\[ \lambda^A_1 = 0.2190, \quad \lambda^A_2 = -0.2190 \]

\[ \lambda^B_1 = -3.3756, \quad \lambda^B_2 = 3.3756 \]

\[ \lambda^{AB}_{11} = 0.1114, \quad \lambda^{AB}_{12} = -0.1114, \quad \lambda^{AB}_{21} = -0.1114, \quad \lambda^{AB}_{22} = 0.1114 \]
## 15.6 Hitung Odds Ratio dan Interval Kepercayaan

\[ \text{OR} = \frac{n_{11}n_{22}}{n_{12}n_{21}} = \frac{28 \times 658}{656 \times 18} = \frac{18424}{11808} = 1.5609 \]

Log odds ratio:

\[ \log(\text{OR}) = \log(1.5609) = 0.4455 \]

Standard error (SE):

\[ SE = \sqrt{\frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}}} \]

\[ SE = \sqrt{ \frac{1}{28} + \frac{1}{656} + \frac{1}{18} + \frac{1}{658} } = \sqrt{0.0357 + 0.0015 + 0.0556 + 0.0015} = \sqrt{0.0943} \] \[ = 0.3071 \]

95% Confidence Interval for log(OR):

\[ \log(OR) \pm 1.96 \times SE = 0.4455 \pm 1.96 \times 0.3071 \]

\[ = (0.4455 - 0.6019,\; 0.4455 + 0.6019) = (-0.1564,\; 1.0474) \]

Back-transform to get CI for OR:

\[ \text{Lower} = \exp(-0.1564) = 0.8555 \]

\[ \text{Upper} = \exp(1.0474) = 2.8501 \]

Jadi, OR = 1.56 (95% CI: 0.86 – 2.85)

Fitting Model Log Linear dengan R

# Data 2x2
tabel <- matrix(c(28, 656, 18, 658), nrow = 2, byrow = TRUE)
colnames(tabel) <- c("Infeksi", "Tidak")
rownames(tabel) <- c("Placebo", "Aspirin")
tabel
##         Infeksi Tidak
## Placebo      28   656
## Aspirin      18   658
# Konversi jadi data.frame untuk model GLM
data <- as.data.frame(as.table(tabel))
colnames(data) <- c("Obat", "Infeksi", "Freq")
data
# Fit model tanpa interaksi
fit_no_inter <- glm(Freq ~ Obat + Infeksi, family = poisson(), data = data)
summary(fit_no_inter)
## 
## Call:
## glm(formula = Freq ~ Obat + Infeksi, family = poisson(), data = data)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   3.14136    0.14989  20.958   <2e-16 ***
## ObatAspirin  -0.01176    0.05423  -0.217    0.828    
## InfeksiTidak  3.35219    0.15000  22.348   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 1485.5614  on 3  degrees of freedom
## Residual deviance:    2.1474  on 1  degrees of freedom
## AIC: 34.713
## 
## Number of Fisher Scoring iterations: 4
# Fit model dengan interaksi
fit_inter <- glm(Freq ~ Obat * Infeksi, family = poisson(), data = data)
summary(fit_inter)
## 
## Call:
## glm(formula = Freq ~ Obat * Infeksi, family = poisson(), data = data)
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                3.3322     0.1890  17.632   <2e-16 ***
## ObatAspirin               -0.4418     0.3021  -1.462    0.144    
## InfeksiTidak               3.1540     0.1930  16.344   <2e-16 ***
## ObatAspirin:InfeksiTidak   0.4449     0.3071   1.449    0.147    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 1.4856e+03  on 3  degrees of freedom
## Residual deviance: 1.3323e-15  on 0  degrees of freedom
## AIC: 34.565
## 
## Number of Fisher Scoring iterations: 3

Interpretasi Parameter

  • Parameter utama (intercept) menunjukkan rata-rata log frekuensi sel referensi (Placebo + Infeksi).

  • Efek “Obat” dan “Infeksi” menunjukkan perbedaan log-frekuensi antara kategori. Efek infeksi signifikan; efek obat tidak.

  • Efek interaksi tidak signifikan, namun positif, menunjukkan bahwa efek Aspirin terhadap infeksi mungkin berbeda tergantung status infeksi, meski bukti statistiknya belum kuat

Analisis Data Tabel Kontingensi 2x3

Sebuah studi meneliti hubungan antara kejadian Myocardial Infarction (MI) dan tingkat merokok (Smoking Level):

0 batang 1–24 batang >25 batang
kontrol 25 25 12
Infeksi Myocardial 0 1 3

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^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]

dengan:

  • \(\mu_{ij}\): ekspektasi frekuensi untuk baris ke-\(i\), kolom ke-\(j\)
  • \(A\): Grup (Kontrol atau IM), \(i = 1, 2\)
  • \(B\): Kategori Merokok (0, 1–24, >25 batang), \(j = 1, 2, 3\)

Constraint: \[ \sum_i \lambda^A_i = 0,\quad \sum_j \lambda^B_j = 0,\quad \sum_i \lambda^{AB}_{ij} = 0,\quad \sum_j \lambda^{AB}_{ij} = 0 \]

Secara eksplisit:

\[ \log(\mu_{ij}) = \lambda + \begin{cases} \lambda^A_1 & \text{(Kontrol)} \\ \lambda^A_2 & \text{(Infeksi Myocardial)} \end{cases} + \begin{cases} \lambda^B_1 & \text{(Merokok = 0)} \\ \lambda^B_2 & \text{(Merokok = 1–24)} \\ \lambda^B_3 & \text{(Merokok > 25)} \end{cases} + \lambda^{AB}_{ij} \text{ (jika ada interaksi)} \]

Fitting Model Log Linear di R

tabel_mi <- matrix(c(25, 25, 12,
                     0,  1,  3), 
                   nrow = 2, byrow = TRUE)
colnames(tabel_mi) <- c("Merokok0", "Merokok1_24", "Merokok_>25")
rownames(tabel_mi) <- c("Kontrol", "IM")
tabel_mi
##         Merokok0 Merokok1_24 Merokok_>25
## Kontrol       25          25          12
## IM             0           1           3
# Ubah menjadi data.frame untuk glm
data_mi <- as.data.frame(as.table(tabel_mi))
colnames(data_mi) <- c("Grup", "Merokok", "Frekuensi")
data_mi
# Model log-linear tanpa interaksi (asumsi independen)
model_tanpa_interaksi <- glm(Frekuensi ~ Grup + Merokok, 
                             family = poisson(), 
                             data = data_mi)
summary(model_tanpa_interaksi)
## 
## Call:
## glm(formula = Frekuensi ~ Grup + Merokok, family = poisson(), 
##     data = data_mi)
## 
## Coefficients:
##                    Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         3.15636    0.20243  15.592  < 2e-16 ***
## GrupIM             -2.74084    0.51588  -5.313 1.08e-07 ***
## MerokokMerokok1_24  0.03922    0.28011   0.140    0.889    
## MerokokMerokok_>25 -0.51083    0.32660  -1.564    0.118    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 71.5948  on 5  degrees of freedom
## Residual deviance:  6.6901  on 2  degrees of freedom
## AIC: 34.145
## 
## Number of Fisher Scoring iterations: 6
# Model log-linear dengan interaksi (untuk cek asosiasi)
model_dengan_interaksi <- glm(Frekuensi ~ Grup * Merokok, 
                              family = poisson(), 
                              data = data_mi)
summary(model_dengan_interaksi)
## 
## Call:
## glm(formula = Frekuensi ~ Grup * Merokok, family = poisson(), 
##     data = data_mi)
## 
## Coefficients:
##                             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                3.219e+00  2.000e-01  16.094   <2e-16 ***
## GrupIM                    -2.552e+01  4.225e+04  -0.001   0.9995    
## MerokokMerokok1_24         1.402e-15  2.828e-01   0.000   1.0000    
## MerokokMerokok_>25        -7.340e-01  3.512e-01  -2.090   0.0366 *  
## GrupIM:MerokokMerokok1_24  2.230e+01  4.225e+04   0.001   0.9996    
## GrupIM:MerokokMerokok_>25  2.414e+01  4.225e+04   0.001   0.9995    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 7.1595e+01  on 5  degrees of freedom
## Residual deviance: 4.1224e-10  on 0  degrees of freedom
## AIC: 31.455
## 
## Number of Fisher Scoring iterations: 20

Interpretasi

Karena interaksi tidak signifikan, maka distribusi kebiasaan merokok serupa antara kelompok MI dan kontrol, atau dengan kata lain tidak ada bukti asosiasi kuat antara MI dan merokok dalam kategori ini.

Model Log Linear 3 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).

2. 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.

3. 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.

4. Model Joint Independence

  • Independensi antara X & Y:

\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]

  • Independensi antara X & Z:

\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \]

  • Independensi antara Y & Z:

\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \]

5. 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 Tiga 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

Soal

Tabel berikut menyajikan data dari survei mengenai penggunaan alkohol, rokok, dan ganja oleh siswa SMA. 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.

Tabel Data Survei

Alcohol Use Cigarette Use Marijuana Use: Yes Marijuana Use: No
Yes Yes 911 538
Yes No 44 456
No Yes 3 43
No No 2 279

Keterangan:

  • Penggunaan Alkohol: Ya, Tidak
  • Penggunaan Rokok: Ya, Tidak
  • Penggunaan Ganja: Ya, Tidak

Analisis Log-Linear untuk Tabel Tiga Arah

Package yang Digunakan

library("epitools")
library("DescTools")
library("lawstat")
## Warning: package 'lawstat' was built under R version 4.4.3

Input Data

# Input data sesuai tabel penggunaan Alkohol, Rokok, dan Ganja
alcohol   <- factor(rep(c("1Yes", "2No"), each = 4))          # 2 levels × 4 = 8
cigarette <- factor(rep(c("1Yes", "1Yes", "2No", "2No"), 2))  # Ulangin karena nested
marijuana <- factor(rep(c("1Yes", "2No"), times = 4))         # Yes, No → berulang

counts <- c(911, 538, 44, 456, 3, 43, 2, 279)

data <- data.frame(
  Alcohol   = alcohol,
  Cigarette = cigarette,
  Marijuana = marijuana,
  Frekuensi = counts
)

data

Membentuk Tabel Kontigensi 3 Arah

# Buat table 3 dimensi dari data alkohol, rokok, ganja
table3d <- xtabs(Frekuensi ~ Alcohol + Cigarette + Marijuana, data = data)

# Tampilkan dengan format flat table
ftable(table3d)
##                   Marijuana 1Yes 2No
## Alcohol Cigarette                   
## 1Yes    1Yes                 911 538
##         2No                   44 456
## 2No     1Yes                   3  43
##         2No                    2 279

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

data$Alcohol   <- relevel(data$Alcohol, ref = "2No")
data$Cigarette <- relevel(data$Cigarette, ref = "2No")
data$Marijuana <- relevel(data$Marijuana, ref = "2No")

Model Saturated Model log-linear 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 ~ Alcohol + Cigarette + Marijuana +
             Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana +
               Alcohol:Cigarette:Marijuana,
             family = poisson(link = "log"), data = data)

summary(model_saturated)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Alcohol:Marijuana + Cigarette:Marijuana + Alcohol:Cigarette:Marijuana, 
##     family = poisson(link = "log"), data = data)
## 
## Coefficients:
##                                         Estimate Std. Error z value Pr(>|z|)
## (Intercept)                              5.63121    0.05987  94.060  < 2e-16
## Alcohol1Yes                              0.49128    0.07601   6.464 1.02e-10
## Cigarette1Yes                           -1.87001    0.16383 -11.414  < 2e-16
## Marijuana1Yes                           -4.93806    0.70964  -6.959 3.44e-12
## Alcohol1Yes:Cigarette1Yes                2.03538    0.17576  11.580  < 2e-16
## Alcohol1Yes:Marijuana1Yes                2.59976    0.72698   3.576 0.000349
## Cigarette1Yes:Marijuana1Yes              2.27548    0.92746   2.453 0.014149
## Alcohol1Yes:Cigarette1Yes:Marijuana1Yes  0.58951    0.94236   0.626 0.531600
##                                            
## (Intercept)                             ***
## Alcohol1Yes                             ***
## Cigarette1Yes                           ***
## Marijuana1Yes                           ***
## Alcohol1Yes:Cigarette1Yes               ***
## Alcohol1Yes:Marijuana1Yes               ***
## Cigarette1Yes:Marijuana1Yes             *  
## Alcohol1Yes:Cigarette1Yes:Marijuana1Yes    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2.8515e+03  on 7  degrees of freedom
## Residual deviance: 1.3056e-13  on 0  degrees of freedom
## AIC: 65.043
## 
## Number of Fisher Scoring iterations: 3
exp(model_saturated$coefficients)
##                             (Intercept)                             Alcohol1Yes 
##                            2.790000e+02                            1.634409e+00 
##                           Cigarette1Yes                           Marijuana1Yes 
##                            1.541219e-01                            7.168459e-03 
##               Alcohol1Yes:Cigarette1Yes               Alcohol1Yes:Marijuana1Yes 
##                            7.655141e+00                            1.346053e+01 
##             Cigarette1Yes:Marijuana1Yes Alcohol1Yes:Cigarette1Yes:Marijuana1Yes 
##                            9.732558e+00                            1.803106e+00

Ringkasan Model

Model ini adalah model log-linear saturated dengan semua efek utama, interaksi dua arah, dan interaksi tiga arah. Model ini memodelkan hubungan antara Penggunaan Alkohol, Penggunaan Rokok, dan Penggunaan Ganja terhadap frekuensi responden.

Hasil Estimasi Koefisien

library(knitr)

signif_stars <- function(pval) {
  if (pval < 0.001) return("***")
  else if (pval < 0.01) return("**")
  else if (pval < 0.05) return("*")
  else if (pval < 0.1) return(".")
  else return("")
}

coef_saturated <- round(summary(model_saturated)$coefficients, 3)

colnames(coef_saturated) <- c("Est", "SE", "z", "p")

Signif <- sapply(coef_saturated[, "p"], signif_stars)
coef_table <- cbind(coef_saturated, Sig = Signif)

kable(coef_table, format = "pipe", caption = "Tabel Koefisien Model Saturated")
Tabel Koefisien Model Saturated
Est SE z p Sig
(Intercept) 5.631 0.06 94.06 0 ***
Alcohol1Yes 0.491 0.076 6.464 0 ***
Cigarette1Yes -1.87 0.164 -11.414 0 ***
Marijuana1Yes -4.938 0.71 -6.959 0 ***
Alcohol1Yes:Cigarette1Yes 2.035 0.176 11.58 0 ***
Alcohol1Yes:Marijuana1Yes 2.6 0.727 3.576 0 ***
Cigarette1Yes:Marijuana1Yes 2.275 0.927 2.453 0.014 *
Alcohol1Yes:Cigarette1Yes:Marijuana1Yes 0.59 0.942 0.626 0.532

Interpretasi Koefisien

(Intercept): Rata-rata log jumlah kasus untuk kategori referensi (tidak minum alkohol, tidak merokok, tidak pakai ganja) adalah 5.631, atau rata-rata kasus sekitar exp(5.631) ≈ 279.

Alcohol1Yes: Mereka yang mengonsumsi alkohol memiliki jumlah kejadian yang 1.63 kali lebih tinggi dibanding yang tidak (exp(0.491) ≈ 1.63), signifikan (p < 0.001).

Cigarette1Yes: Mereka yang merokok memiliki jumlah kejadian sekitar 0.15 kali lebih rendah dari yang tidak merokok (exp(-1.870) ≈ 0.15), signifikan (p < 0.001).

Marijuana1Yes: Mereka yang menggunakan ganja memiliki jumlah kejadian sekitar 0.007 kali lebih rendah dibanding yang tidak (exp(-4.938) ≈ 0.007), signifikan (p < 0.001).

Goodness of Fit

Residual deviance ≈ 0 menandakan model saturated benar-benar fit terhadap data (seluruh variasi data dijelaskan oleh model).

Kesimpulan

  • Konsumsi alkohol, konsumsi rokok, dan konsumsi ganja masing-masing secara signifikan meningkatkan jumlah kejadian.

  • Interaksi antara alkohol dan rokok, antara alkohol dan ganja, serta rokok dan ganja secara signifikan meningkatkan jumlah kejadian.

  • Interaksi tiga zat (alkohol, rokok, dan ganja) tidak menunjukkan pengaruh yang signifikan.

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_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

# Model homogenous (tanpa interaksi tiga arah)
model_homogenous <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana,
              family = poisson(link = "log"), data = data)

summary(model_homogenous)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Alcohol:Marijuana + Cigarette:Marijuana, family = poisson(link = "log"), 
##     data = data)
## 
## Coefficients:
##                             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  5.63342    0.05970  94.361  < 2e-16 ***
## Alcohol1Yes                  0.48772    0.07577   6.437 1.22e-10 ***
## Cigarette1Yes               -1.88667    0.16270 -11.596  < 2e-16 ***
## Marijuana1Yes               -5.30904    0.47520 -11.172  < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes    2.05453    0.17406  11.803  < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes    2.98601    0.46468   6.426 1.31e-10 ***
## Cigarette1Yes:Marijuana1Yes  2.84789    0.16384  17.382  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.46098  on 7  degrees of freedom
## Residual deviance:    0.37399  on 1  degrees of freedom
## AIC: 63.417
## 
## Number of Fisher Scoring iterations: 4

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

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

Langkah-Langkah Pengujian

Hipotesis

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

H1: Ada interaksi tiga arah (model saturated diperlukan)

Hitung Selisih Deviance

# Deviance antar model
Deviance.model <- model_homogenous$deviance - model_saturated$deviance
Deviance.model
## [1] 0.3739859

Hitung Derajat Bebas

# Derajat bebas = db model homogenous - db model saturated
derajat.bebas <- (model_homogenous$df.residual - model_saturated$df.residual)
derajat.bebas
## [1] 1

Chi Square Tabel

chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel
## [1] 3.841459

Keputusan Uji

Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan
## [1] "Terima H0"

Interpretasi Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara Penggunaan Alkohol, Penggunaan Rokok, dan Penggunaan Ganja. Atau dengan kata lain, model homogenous saja sudah cukup.

Rangkuman

Pengujian Ada Tidaknya Interaksi Tiga Arah (Saturated Model vs Homogenous Model)

Hipotesis:

  • \(H_0\): \(\lambda_{ijk}^{XYZ} = 0\) (Tidak ada interaksi tiga arah; model yang terbentuk adalah model homogenous)
  • \(H_1\): \(\lambda_{ijk}^{XYZ} \neq 0\) (Ada interaksi tiga arah; model yang terbentuk adalah model saturated)

Tingkat Signifikansi:

  • \(\alpha = 5\%\)

Statistik Uji:

\[ \Delta \text{Deviance} = \text{Deviance model homogenous} - \text{Deviance model saturated} \] \[ = 0.374 - 0.0 = 0.374 \]

\[ df = df\ \text{model homogenous} - df\ \text{model saturated} = 1 - 0 = 1 \]

Daerah Penolakan:

  • Tolak \(H_0\) jika \(\Delta \text{Deviance} > \chi^2_{0.05, 1} = 3.841\)

Keputusan:

  • Karena \(0.374 < 3.841\), maka terima \(H_0\)

Interpretasi:

  • Pada taraf nyata 5%, belum cukup bukti untuk menolak \(H_0\), atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara Penggunaan Alkohol, Penggunaan Rokok, dan Penggunaan Ganja.

Uji Model Interaksi Dua Arah (Homogenous VS Conditional On X)

Model 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_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} \]

# Model Conditional Association on X
model_conditional_X <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Cigarette + Alcohol:Marijuana,
              family = poisson(link = "log"), data = data)

summary(model_conditional_X)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Alcohol:Marijuana, family = poisson(link = "log"), data = data)
## 
## Coefficients:
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                5.62295    0.06005  93.635   <2e-16 ***
## Alcohol1Yes               -0.08167    0.07810  -1.046    0.296    
## Cigarette1Yes             -1.80971    0.15905 -11.378   <2e-16 ***
## Marijuana1Yes             -4.16511    0.45067  -9.242   <2e-16 ***
## Alcohol1Yes:Cigarette1Yes  2.87373    0.16730  17.178   <2e-16 ***
## Alcohol1Yes:Marijuana1Yes  4.12509    0.45294   9.107   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.46  on 7  degrees of freedom
## Residual deviance:  497.37  on 2  degrees of freedom
## AIC: 558.41
## 
## Number of Fisher Scoring iterations: 5

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

• Hipotesis \[ H_0 : \lambda_{jk}^{YZ} = 0 \] Tidak ada interaksi antara Penggunaan Rokok (y) dan Penggunaan Ganja (z)

\[ H_1 : \lambda_{jk}^{YZ} \ne 0 \] Ada interaksi antara Penggunaan Rokok (y) dan Penggunaan Ganja (z)

• Taraf Signifikansi \[ \alpha = 5\% \]

• Statistik Uji \[ \Delta \text{Deviance} = \text{Deviance model conditional on } X - \text{Deviance model homogenous} \]

\[ = 497.369 - 0.374 = 496.995 \]

\[ db = \text{df model conditional on } X - \text{df model homogenous} \]

\[ = 2 - 1 = 1 \]

• Daerah Penolakan \[ \text{Tolak } H_0 \text{ jika } \Delta \text{Deviance} > \chi^2_{0.05, 1} = 3.841 \] • Keputusan

Karena 496.95 > 3.841, maka tolak \(H_0\)

• Kesimpulan Dengan taraf nyata 5%, ada interaksi yang signifikan antara Penggunaan Rokok (Y) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{jk}^{YZ}\).

Pengujian Selisih Deviance (Conditional on X vs Homogenous)

# Selisih deviance antar model
Deviance.model <- model_conditional_X$deviance - model_homogenous$deviance
Deviance.model
## [1] 496.9953

Hitung Derajat Bebas

# Chi Square tabel dengan alpha = 0.05
derajat.bebas <- (2 - 1)
derajat.bebas
## [1] 1

Nilai Chi-Square Tabel

chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 3.841459

Keputusan Uji

Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Tolak"

Interpretasi

Karena nilai Deviance.model = 996.95 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 bukti yang cukup untuk menolak \(H_0\). Artinya, ada interaksi yang signifkan antara Penggunaan Rokok (Y) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model yang menyertakan parameter interaksi \(\lambda_{jk}^{YZ}\).

Uji Model Interaksi Dua Arah (Homogenous VS Conditional On Y)

Model 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_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{jk}^{YZ} \]

# Model Conditional Association on Y (Cigarette)
model_conditional_Y <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Cigarette + Cigarette:Marijuana,
              family = poisson(link = "log"), data = data)

summary(model_conditional_Y)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Cigarette:Marijuana, family = poisson(link = "log"), data = data)
## 
## Coefficients:
##                             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  5.57765    0.06032  92.463  < 2e-16 ***
## Alcohol1Yes                  0.57625    0.07456   7.729 1.08e-14 ***
## Cigarette1Yes               -2.69414    0.16257 -16.572  < 2e-16 ***
## Marijuana1Yes               -2.77123    0.15199 -18.233  < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes    2.87373    0.16730  17.178  < 2e-16 ***
## Cigarette1Yes:Marijuana1Yes  3.22431    0.16098  20.029  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.461  on 7  degrees of freedom
## Residual deviance:   92.018  on 2  degrees of freedom
## AIC: 153.06
## 
## Number of Fisher Scoring iterations: 6

Pengujian Ada Tidaknya Interaksi antara X dan Z (Homogenous model vs Conditional Association on Y)

• Hipotesis

\[ H_0 : \lambda_{ik}^{XZ} = 0 \] Tidak ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Ganja (z)

\[ H_1 : \lambda_{ik}^{XZ} \ne 0 \] Ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Ganja (z)

• Tingkat Signifikansi

\[ \alpha = 5\% \]

• Statistik Uji

\[ \Delta \text{Deviance} = \text{Deviance model conditional on } Y - \text{Deviance model homogenous} \]

\[ = 92.018 - 0.372 = 91.644 \]

\[ db = \text{df model conditional on } Y - \text{df model homogenous} \]

\[ = 2 - 1 = 1 \]

• Daerah Penolakan

Tolak \(H_0\) jika \(\Delta\text{Deviance} > \chi^2_{0.05, 1} = 3.841\)

• Keputusan

\[ \text{Karena } 91.644 > 3.841, \text{ maka tolak } H_0 \] • Kesimpulan

Dengan taraf nyata 5%, ada interaksi yang signifikan antara Penggunaan Alkohol (X) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{ik}^{XZ}\).

Pengujian Hipotesis Interaksi X dan Z (Conditional on Y vs Homogenous)

# Deviance of Model
Deviance.model <- model_conditional_Y$deviance - model_homogenous$deviance 
Deviance.model
## [1] 91.64437

Hitung Derajat Bebas

derajat.bebas <- (2 - 1)
derajat.bebas
## [1] 1

Nilai Chi-Square Tabel

chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 3.841459

Keputusan Uji

Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Tolak"

Interpretasi

Karena nilai Deviance.model = 91.644 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 bukti yang cukup untuk menolak \(H_0\). Artinya, ada interaksi yang signifkan antara Penggunaan Alkohol (X)) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model yang menyertakan parameter interaksi \(\lambda_{ik}^{XZ}\).

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

Model 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_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

# Model Conditional Association on Z (Marijuana)
model_conditional_Z <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Marijuana + Cigarette:Marijuana,
              family = poisson(link = "log"), data = data)

summary(model_conditional_Z)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Marijuana + 
##     Cigarette:Marijuana, family = poisson(link = "log"), data = data)
## 
## Coefficients:
##                             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  5.19207    0.06088  85.285  < 2e-16 ***
## Alcohol1Yes                  1.12719    0.06412  17.579  < 2e-16 ***
## Cigarette1Yes               -0.23512    0.05551  -4.235 2.28e-05 ***
## Marijuana1Yes               -6.62092    0.47370 -13.977  < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes    4.12509    0.45294   9.107  < 2e-16 ***
## Cigarette1Yes:Marijuana1Yes  3.22431    0.16098  20.029  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.46  on 7  degrees of freedom
## Residual deviance:  187.75  on 2  degrees of freedom
## AIC: 248.8
## 
## Number of Fisher Scoring iterations: 5

Pengujian Ada Tidaknya Interaksi antara X dan Y (Homogenous model vs Conditional Association on Z)

• Hipotesis

\[ H_0 : \lambda_{ij}^{XY} = 0 \] Tidak ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Rokok (y)

\[ H_1 : \lambda_{ij}^{XY} \ne 0 \] Ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Rokok (y)

• Tingkat Signifikansi

\[ \alpha = 5\% \]

• Statistik Uji

\[ \Delta\text{Deviance} = \text{Deviance model conditional on } Z - \text{Deviance model homogenous} \]

\[ = 187.75 - 0.374 = 187.380 \]

\[ db = \text{df model conditional on } Z - \text{df model homogenous} \]

\[ = 2 - 1 = 1 \]

• Daerah Penolakan

Tolak \(H_0\) jika \(\Delta\text{Deviance} > \chi^2_{0.05, 1} = 3.841\)

• Keputusan

\[ \text{Karena } 187.380 > 3.841, \text{ maka tolak } H_0 \] • Kesimpulan

Dengan taraf nyata 5%, Ada interaksi yang signifikan antara Penggunaan Alkohol (X) dan Penggunaan Rokok (Y)). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{ij}^{XY}\).

Pengujian Hipotesis Interaksi X dan Y (Conditional on Z vs Homogenous)

# Deviance of Model
Deviance.model <- model_conditional_Z$deviance - model_homogenous$deviance 
Deviance.model
## [1] 187.3803

Hitung Derajat Bebas

derajat.bebas <- (2 - 1)
derajat.bebas
## [1] 1

Nilai Chi-Square Tabel

chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 3.841459

Keputusan Uji

Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Tolak"

Interpretasi

Karena nilai Deviance.model = 187.376 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 bukti yang cukup untuk menolak \(H_0\). Artinya, ada interaksi yang signifkan antara Penggunaan Alkohol (X)) dan Penggunaan Rokok (Y)). Dengan kata lain, model yang terbentuk adalah model yang menyertakan parameter interaksi \(\lambda_{ij}^{XY}\).

Pemilihan Model Terbaik

Ringkasan Model Log Linier

Model log linier dengan berbagai kombinasi parameter:

library(knitr)

model_table <- data.frame(
  Model = c(
    "Saturated",
    "Homogenous",
    "Conditional on X",
    "Conditional on Y",
    "Conditional on Z"
  ),
  Parameter = c(
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ", 
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ",             
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ",                        
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ",                        
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ"                        
  ),
  Deviance = c(0.000, 0.374, 497.37, 92.018, 187.75),
  Param = c(8, 7, 6, 6, 6),
  df = c(0, 1, 1, 1, 1),
  AIC = c(65.043, 63.417, 558.41, 153.06, 248.80)
)

kable(
  model_table,
  format = "pipe",
  caption = "Ringkasan Model Log-Linear"
)
Ringkasan Model Log-Linear
Model Parameter Deviance Param df AIC
Saturated λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ 0.000 8 0 65.043
Homogenous λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ 0.374 7 1 63.417
Conditional on X λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ 497.370 6 1 558.410
Conditional on Y λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ 92.018 6 1 153.060
Conditional on Z λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ 187.750 6 1 248.800

Ringkasan Pengujian Interaksi 3 Arah dan 2 Arah

library(knitr)
uji_interaksi <- data.frame(
  Interaksi = c("XYZ", "YZ", "XZ", "XY"),
  Pengujian = c(
    "Saturated vs Homogenous",
    "Conditional on X vs Homogenous",
    "Conditional on Y vs Homogenous",
    "Conditional on Z vs Homogenous"
  ),
  "Δ Deviance" = c(0.374, 496.995, 91.644, 187.380),
  df = c(2, 1, 1, 1),
  ChiSqTabel = c(5.991, 3.841, 3.841, 3.841),
  Ket = c("Tidak ada interaksi", "Ada interaksi", "Ada interaksi", "Ada interaksi"),
  Keputusan = c("Terima $H_0$", "Tolak $H_0$", "Tolak $H_0$", "Tolak $H_0$")
)

# Tampilkan dengan kable
kable(uji_interaksi, format = "pipe", caption = "Ringkasan Uji Interaksi (Chi-Square)")
Ringkasan Uji Interaksi (Chi-Square)
Interaksi Pengujian Δ.Deviance df ChiSqTabel Ket Keputusan
XYZ Saturated vs Homogenous 0.374 2 5.991 Tidak ada interaksi Terima \(H_0\)
YZ Conditional on X vs Homogenous 496.995 1 3.841 Ada interaksi Tolak \(H_0\)
XZ Conditional on Y vs Homogenous 91.644 1 3.841 Ada interaksi Tolak \(H_0\)
XY Conditional on Z vs Homogenous 187.380 1 3.841 Ada interaksi Tolak \(H_0\)

Kesimpulan Pemilihan Model Terbaik

Dari hasil di atas diketahui bahwa asosiasi yang nyata terdapat pada model Homogenous dan ketiga model Conditional. Karena itulah, model akan dipilih berdasarkan nilai AIC yang paling kecil, dimana AIC terkecil terdapat pada model Homogenous yaitu sebesar 63.417. Oleh karena itu, model terbaik yaitu:

\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

Model terbaik adalah model log-linear tanpa interaksi tiga arah dan memuat semua interaksi dua arah.

Model Terbaik

Model terbaik dipilih berdasarkan nilai AIC paling kecil, yaitu model Homogenous.

\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

# Model Terbaik
best_model <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana,
              family = poisson(link = "log"), data = data)

summary(best_model)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Alcohol:Marijuana + Cigarette:Marijuana, family = poisson(link = "log"), 
##     data = data)
## 
## Coefficients:
##                             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  5.63342    0.05970  94.361  < 2e-16 ***
## Alcohol1Yes                  0.48772    0.07577   6.437 1.22e-10 ***
## Cigarette1Yes               -1.88667    0.16270 -11.596  < 2e-16 ***
## Marijuana1Yes               -5.30904    0.47520 -11.172  < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes    2.05453    0.17406  11.803  < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes    2.98601    0.46468   6.426 1.31e-10 ***
## Cigarette1Yes:Marijuana1Yes  2.84789    0.16384  17.382  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.46098  on 7  degrees of freedom
## Residual deviance:    0.37399  on 1  degrees of freedom
## AIC: 63.417
## 
## Number of Fisher Scoring iterations: 4

Dari summary model diatas terlihat bahwa best model memiliki AIC yang lebih rendah dibandingkan saturated dan conditional model.

Interpretasi Koefisien Model Terbaik

# Interpretasi koefisien model terbaik
data.frame(
koef = best_model$coefficients,
exp_koef = exp(best_model$coefficients)
)

Interpretasi Koefisien Model Terbaik

\[ \exp(\lambda^{X}_{1A}) = \exp(0.488) = 1.629 \Rightarrow \text{nilai odds} \]

Tanpa memperhatikan penggunaan rokok dan ganja, mereka yang mengonsumsi alkohol memiliki jumlah kasus 1.63 kali lebih besar dibanding yang tidak mengonsumsi alkohol

\[ \exp(\lambda^{Y}_{1C}) = \exp(-1.887) = 0.152 \Rightarrow \text{nilai odds} \]

Tanpa memperhatikan penggunaan alkohol dan ganja, mereka yang merokok memiliki jumlah kasus 0.15 kali lebih kecil dibanding yang tidak merokok

\[ \exp(\lambda^{Z}_{1M}) = \exp(-5.309) = 0.004 \Rightarrow \text{nilai odds} \]

Tanpa memperhatikan penggunaan alkohol dan rokok, mereka yang menggunakan ganja memiliki jumlah kasus hanya sekitar 0.004 kali dari yang tidak menggunakan ganja

\[ \exp(\lambda^{Z}_{1A,1C}) = \exp(2.055) = 7.807 \Rightarrow \text{nilai odds} \]

Responden yang menggunakan Alkohol dan Rokok memiliki jumlah kasus 7.8 kali lebih besar dibanding jika hanya salah satu atau tidak sama sekali.

\[ \exp(\lambda^{XY}_{1A,1M}) = \exp(2.986) = 19.806 \Rightarrow \text{nilai odds ratio} \]

Responden yang menggunakan Alkohol dan Ganja memiliki jumlah kasus 19.8 kali lebih besar dibanding jika hanya salah satu atau tidak sama sekali.

\[ \exp(\lambda^{XY}_{1C, 1M}) = \exp(2.848) = 17.253 \Rightarrow \text{nilai odds ratio} \] Responden yang menggunakan Rokok dan Ganja memiliki jumlah kasus 17.25 kali lebih besar dibanding jika hanya salah satu atau tidak sama sekali.

Nilai Dugaan Model Terbaik

data.frame(
  Alcohol   = data$Alcohol,
  Cigarette = data$Cigarette,
  Marijuana = data$Marijuana,
  Frekuensi = data$Frekuensi,
  fitted    = fitted(model_homogenous)
)

Perhitungan Manual Nilai Dugaan (Fitted Value) Model Terbaik

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

\[ \hat{\mu}_{111} = \exp(5.631 + 0.491 - 1.870 - 4.938 + 2.035 + 2.600 + 2.275 + 0.590) = \]

\[ \exp(6.811) = 910.383 \]

\[ \hat{\mu}_{112} = \exp(5.631 + 0.491 - 1.870 + 0 + 2.035 + 0 + 0 + 0) = \]

\[ \exp(6.287) = 538.617 \]

\[ \hat{\mu}_{121} = \exp(5.631 + 0.491 + 0 - 4.938 + 0 + 2.600 + 0 + 0) = \]

\[ \exp(3.799) = 44.617 \]

\[ \hat{\mu}_{122} = \exp(5.631 + 0.491 + 0 + 0 + 0 + 0 + 0 + 0) = \]

\[ \exp(6.122) = 455.383 \]

\[ \hat{\mu}_{211} = \exp(5.631 + 0 + -1.870 - 4.938 + 0 + 0 + 2.275 + 0) = \]

\[ \exp(1.285) = 3.617 \]

\[ \hat{\mu}_{212} = \exp(5.631 + 0 - 1.870 + 0 + 0 + 0 + 0 + 0) = \]

\[ \exp(3.746) = 42.383 \]

\[ \hat{\mu}_{221} = \exp(5.631 + 0 + 0 - 4.938 + 0 + 0 + 0 + 0) = \]

\[ \exp(1.285) = 1.383 \]

\[ \hat{\mu}_{222} = \exp(5.631 + 0 + 0 + 0 + 0 + 0 + 0 + 0) = \]

\[ \exp(5.631) = 279.617 \]

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

Interpretasi:

Nilai fitted values merupakan prediksi dari model yang sudah mempertimbangkan interaksi dua arah. Nilai ini diperoleh dari estimasi model log-linear yang menyertakan semua efek utama dan semua interaksi dua arah. Dari perbandingan antara prediksi model dan nilai asli, didapatkan bahwa nilai tersebut saling mirip serta tidak terdapat perbedaan yang signifikan. Oleh karena itu, dapat disimpulkan bahwa model dapat memprediksi data dengan sangat baik.

Studi Kasus 1

Soal

Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.

Penggunaan Sabuk Pengaman Cedera Fatal Cedera Tidak Fatal Total
Tidak memakai sabuk pengaman 1,601 162,527 164,128
Memakai sabuk pengaman 510 412,368 412,878
Total 2,111 574,895 577,006

Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988

Hitunglah peluang bersama, peluang marjinal, dan peluang bersyarat beserta uji independensi dan ukuran asosiasi

Jawab

Peluang Langkah 1 : Hitung Peluang Bersama

  • \[P(\text{No Seatbelt, Fatal}) = \frac{1601}{577006} = 0.0028\]
  • \[P(\text{No Seatbelt, NonFatal}) = \frac{162527}{577006} = 0.2816\]
  • \[P(\text{Seatbelt, Fatal}) = \frac{510}{577006} = 0.0009\]
  • \[P(\text{Seatbelt, NonFatal}) = \frac{412368}{577006} = 0.7147\]

Langkah 2 : Hitung Peluang Marginal

  • \[P(\text{No Seatbelt}) = \frac{164128}{577006} = 0.2844\]
  • \[P(\text{Seatbelt}) = \frac{412878}{577006} = 0.7156\]
  • \[P(\text{Fatal}) = \frac{2111}{577006} = 0.0037\]
  • \[P(\text{NonFatal}) = \frac{574895}{577006} = 0.9963\]

Langkah 3 : Hitung Peluang Bersyarat

  • \[ P(\text{Fatal}|\text{No Seatbelt}) = \frac{1601}{164128} = 0.0098 \]

  • \[ P(\text{Fatal}|\text{Seatbelt}) = \frac{510}{412878} = 0.0012 \]

  • \[ P(\text{Seatbelt}|\text{Fatal}) = \frac{510}{2111} = 0.2416 \]

  • \[ P(\text{Seatbelt}|\text{NonFatal}) = \frac{412368}{574895} = 0.7173 \]

# Membuat tabel 2x2
data <- matrix(c(1601, 162527, 510, 412368), nrow = 2, byrow = TRUE)
rownames(data) <- c("Tidak Memakai Sabuk", "Memakai Sabuk")
colnames(data) <- c("Cedera Fatal", "Cedera Tidak Fatal")

# Total keseluruhan
total <- sum(data)

# Peluang bersama (joint probabilities)
joint_prob <- data / total
joint_prob
##                     Cedera Fatal Cedera Tidak Fatal
## Tidak Memakai Sabuk  0.002774668          0.2816730
## Memakai Sabuk        0.000883873          0.7146685
# Peluang marginal (baris dan kolom)
marginal_row <- rowSums(data) / total
marginal_col <- colSums(data) / total
marginal_row
## Tidak Memakai Sabuk       Memakai Sabuk 
##           0.2844476           0.7155524
marginal_col
##       Cedera Fatal Cedera Tidak Fatal 
##        0.003658541        0.996341459
# Peluang bersyarat
# P(Cedera Fatal | Tidak Memakai Sabuk)
p_fatal_given_no_seatbelt <- data[1, 1] / sum(data[1, ])
# P(Cedera Fatal | Memakai Sabuk)
p_fatal_given_seatbelt <- data[2, 1] / sum(data[2, ])
# P(Memakai Sabuk | Cedera Fatal)
p_seatbelt_given_fatal <- data[2, 1] / sum(data[, 1])
# P(Memakai Sabuk | Cedera Tidak Fatal)
p_seatbelt_given_nonfatal <- data[2, 2] / sum(data[, 2])

p_fatal_given_no_seatbelt
## [1] 0.009754582
p_fatal_given_seatbelt
## [1] 0.001235232
p_seatbelt_given_fatal
## [1] 0.2415917
p_seatbelt_given_nonfatal
## [1] 0.7172927

Interpretasi

Peluang Bersama:

  • 0.00277 → Peluang seseorang tidak memakai sabuk dan cedera fatal.

  • 0.00088 → Peluang seseorang memakai sabuk dan cedera fatal.

  • 0.28167 → Peluang seseorang tidak memakai sabuk dan tidak fatal.

  • 0.71467 → Peluang seseorang memakai sabuk dan tidak fatal.

Peluang Marginal:

  • 0.2844 → Peluang seseorang tidak memakai sabuk secara keseluruhan.

  • 0.7156 → Peluang seseorang memakai sabuk secara keseluruhan.

  • 0.00366 → Peluang mengalami cedera fatal secara umum.

  • 0.9963 → Peluang mengalami cedera tidak fatal secara umum.

Peluang Bersyarat:

Dikarenakan \[ \mathit{P}(\text{Fatal} \mid \text{No Seatbelt}) > \mathit{P}(\text{Fatal} \mid \text{Seatbelt}) \] maka dapat disimpulkan bahwa penggunaan sabuk pengaman dapat menurunkan risiko fatalitas dalam kecelakaan.

Uji Independensi

# Membuat data tabel 2x2
data <- matrix(c(1601, 162527, 510, 412368), nrow = 2, byrow = TRUE)

# Menambahkan nama baris dan kolom
rownames(data) <- c("Tidak Memakai Sabuk", "Memakai Sabuk")
colnames(data) <- c("Cedera Fatal", "Cedera Tidak Fatal")

# Menampilkan data tabel
data
##                     Cedera Fatal Cedera Tidak Fatal
## Tidak Memakai Sabuk         1601             162527
## Memakai Sabuk                510             412368
# Menghitung uji 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 = 2336.1, df = 1, p-value < 2.2e-16

Interpretasi:

Dengan taraf signifikansi 5%, didapatkan p-value < 0.05, maka kita menolak H0 (hipotesis nol). Artinya, ada hubungan yang signifikan secara statistik antara penggunaan sabuk pengaman dan tingkat fatalitas kecelakaan

Ukuran Asosiasi

  • Risk Difference \[ RD = \frac{a}{a + b} - \frac{c}{c + d} = \frac{1601}{164128} - \frac{510}{412878} = 0.0098 - 0.0012 = 0.0086 \]

  • Relative Risk \[ \frac{P(\text{Fatal} \mid \text{No Seatbelt})}{P(\text{Fatal} \mid \text{Seatbelt})} = \frac{\frac{1601}{164128}}{\frac{510}{412878}} = \frac{0.009759}{0.001235} \approx 7.897 \]

  • Odds Ratio \[ OR = \frac{a \cdot d}{b \cdot c} = \frac{1601 \cdot 412368}{162527 \cdot 510} = \frac{660201168}{82988770} = 7.96 \]

# Data
a <- 1601
b <- 162527
c <- 510
d <- 412368

# 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.00851935
## 
## $Relative_Risk
## [1] 7.896965
## 
## $Odds_Ratio
## [1] 7.964905

Interpretasi :

  • RD = 0.0086 → Perbedaan risiko kecelakaan fatal sebesar 0.86% lebih tinggi pada pengendara yang tidak pakai sabuk.

  • RR = 7.90 → Pengendara yang tidak memakai sabuk memiliki 8.17 kali risiko lebih besar dalam mengalami kecelakaan fatal.

  • OR = 7.96 → Peluang kecelakaan fatal 7.96 kali lebih besar pada pengendara yang tidak memakai sabuk.

Studi Kasus 2

Soal

Tabel berikut menyajikan data dari survei mengenai penggunaan alkohol, rokok, dan ganja oleh mahasiswa. Tentukan model log-linear terbaik yang akan digunakan dan tunjukkan bagaimana nilai yang diprediksi dari model menggambarkan asosiasi pada data tersebut

Alcohol Use Cigarette Use Marijuana Use: Yes Marijuana Use: No
Yes Yes 911 538
Yes No 44 456
No Yes 3 43
No No 2 279

Sumber: Wright State University School of Medicine & United Health Services. (1992). Survey of substance use among high school seniors in a nonurban area near Dayton, Ohio. Dayton, OH

Keterangan:

  • Penggunaan Alkohol: Ya, Tidak
  • Penggunaan Rokok: Ya, Tidak
  • Penggunaan Ganja: Ya, Tidak

Jawab

Data

library("epitools")
library("DescTools")
library("lawstat")
# Input data sesuai tabel penggunaan Alkohol, Rokok, dan Ganja
alcohol   <- factor(rep(c("1Yes", "2No"), each = 4))          # 2 levels × 4 = 8
cigarette <- factor(rep(c("1Yes", "1Yes", "2No", "2No"), 2))  # Ulangin karena nested
marijuana <- factor(rep(c("1Yes", "2No"), times = 4))         # Yes, No → berulang

counts <- c(911, 538, 44, 456, 3, 43, 2, 279)

data <- data.frame(
  Alcohol   = alcohol,
  Cigarette = cigarette,
  Marijuana = marijuana,
  Frekuensi = counts
)

#Relevel Data
data$Alcohol   <- relevel(data$Alcohol, ref = "2No")
data$Cigarette <- relevel(data$Cigarette, ref = "2No")
data$Marijuana <- relevel(data$Marijuana, ref = "2No")

# Buat table 3 dimensi dari data alkohol, rokok, ganja
table3d <- xtabs(Frekuensi ~ Alcohol + Cigarette + Marijuana, data = data)

# Tampilkan dengan format flat table
ftable(table3d)
##                   Marijuana 2No 1Yes
## Alcohol Cigarette                   
## 2No     2No                 279    2
##         1Yes                 43    3
## 1Yes    2No                 456   44
##         1Yes                538  911

Perbandingan Model

Model Saturated

Model log-linear 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 ~ Alcohol + Cigarette + Marijuana +
             Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana +
               Alcohol:Cigarette:Marijuana,
             family = poisson(link = "log"), data = data)

summary(model_saturated)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Alcohol:Marijuana + Cigarette:Marijuana + Alcohol:Cigarette:Marijuana, 
##     family = poisson(link = "log"), data = data)
## 
## Coefficients:
##                                         Estimate Std. Error z value Pr(>|z|)
## (Intercept)                              5.63121    0.05987  94.060  < 2e-16
## Alcohol1Yes                              0.49128    0.07601   6.464 1.02e-10
## Cigarette1Yes                           -1.87001    0.16383 -11.414  < 2e-16
## Marijuana1Yes                           -4.93806    0.70964  -6.959 3.44e-12
## Alcohol1Yes:Cigarette1Yes                2.03538    0.17576  11.580  < 2e-16
## Alcohol1Yes:Marijuana1Yes                2.59976    0.72698   3.576 0.000349
## Cigarette1Yes:Marijuana1Yes              2.27548    0.92746   2.453 0.014149
## Alcohol1Yes:Cigarette1Yes:Marijuana1Yes  0.58951    0.94236   0.626 0.531600
##                                            
## (Intercept)                             ***
## Alcohol1Yes                             ***
## Cigarette1Yes                           ***
## Marijuana1Yes                           ***
## Alcohol1Yes:Cigarette1Yes               ***
## Alcohol1Yes:Marijuana1Yes               ***
## Cigarette1Yes:Marijuana1Yes             *  
## Alcohol1Yes:Cigarette1Yes:Marijuana1Yes    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2.8515e+03  on 7  degrees of freedom
## Residual deviance: 1.3056e-13  on 0  degrees of freedom
## AIC: 65.043
## 
## Number of Fisher Scoring iterations: 3

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_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

# Model homogenous
model_homogenous <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana,
              family = poisson(link = "log"), data = data)

summary(model_homogenous)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Alcohol:Marijuana + Cigarette:Marijuana, family = poisson(link = "log"), 
##     data = data)
## 
## Coefficients:
##                             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  5.63342    0.05970  94.361  < 2e-16 ***
## Alcohol1Yes                  0.48772    0.07577   6.437 1.22e-10 ***
## Cigarette1Yes               -1.88667    0.16270 -11.596  < 2e-16 ***
## Marijuana1Yes               -5.30904    0.47520 -11.172  < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes    2.05453    0.17406  11.803  < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes    2.98601    0.46468   6.426 1.31e-10 ***
## Cigarette1Yes:Marijuana1Yes  2.84789    0.16384  17.382  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.46098  on 7  degrees of freedom
## Residual deviance:    0.37399  on 1  degrees of freedom
## AIC: 63.417
## 
## Number of Fisher Scoring iterations: 4

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

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

Hipotesis

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

H1: Ada interaksi tiga arah (model saturated diperlukan)

Selisih Deviance

# Deviance antar model
Deviance.model <- model_homogenous$deviance - model_saturated$deviance
Deviance.model
## [1] 0.3739859

Derajat Bebas

# Derajat bebas = db model homogenous - db model saturated
derajat.bebas <- (model_homogenous$df.residual - model_saturated$df.residual)
derajat.bebas
## [1] 1

Chi Square Tabel

chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel
## [1] 3.841459

Keputusan Uji

Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan
## [1] "Terima H0"

Interpretasi

Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara Penggunaan Alkohol, Penggunaan Rokok, dan Penggunaan Ganja. Atau dengan kata lain, model homogenous saja sudah cukup.

Model 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_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} \]

# Model Conditional Association on X
model_conditional_X <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Cigarette + Alcohol:Marijuana,
              family = poisson(link = "log"), data = data)

summary(model_conditional_X)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Alcohol:Marijuana, family = poisson(link = "log"), data = data)
## 
## Coefficients:
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                5.62295    0.06005  93.635   <2e-16 ***
## Alcohol1Yes               -0.08167    0.07810  -1.046    0.296    
## Cigarette1Yes             -1.80971    0.15905 -11.378   <2e-16 ***
## Marijuana1Yes             -4.16511    0.45067  -9.242   <2e-16 ***
## Alcohol1Yes:Cigarette1Yes  2.87373    0.16730  17.178   <2e-16 ***
## Alcohol1Yes:Marijuana1Yes  4.12509    0.45294   9.107   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.46  on 7  degrees of freedom
## Residual deviance:  497.37  on 2  degrees of freedom
## AIC: 558.41
## 
## Number of Fisher Scoring iterations: 5

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

• Hipotesis \[ H_0 : \lambda_{jk}^{YZ} = 0 \] Tidak ada interaksi antara Penggunaan Rokok (y) dan Penggunaan Ganja (z)

\[ H_1 : \lambda_{jk}^{YZ} \ne 0 \] Ada interaksi antara Penggunaan Rokok (y) dan Penggunaan Ganja (z)

• Taraf Signifikansi \[ \alpha = 5\% \]

• Statistik Uji \[ \Delta \text{Deviance} = \text{Deviance model conditional on } X - \text{Deviance model homogenous} \]

\[ = 497.369 - 0.374 = 496.995 \]

\[ db = \text{df model conditional on } X - \text{df model homogenous} \]

\[ = 2 - 1 = 1 \]

• Daerah Penolakan \[ \text{Tolak } H_0 \text{ jika } \Delta \text{Deviance} > \chi^2_{0.05, 1} = 3.841 \] • Keputusan

Karena 496.95 > 3.841, maka tolak \(H_0\)

• Kesimpulan Dengan taraf nyata 5%, ada interaksi yang signifikan antara Penggunaan Rokok (Y) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{jk}^{YZ}\).

Model 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_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{jk}^{YZ} \]

# Model Conditional Association on Y (Cigarette)
model_conditional_Y <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Cigarette + Cigarette:Marijuana,
              family = poisson(link = "log"), data = data)

summary(model_conditional_Y)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette + 
##     Cigarette:Marijuana, family = poisson(link = "log"), data = data)
## 
## Coefficients:
##                             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  5.57765    0.06032  92.463  < 2e-16 ***
## Alcohol1Yes                  0.57625    0.07456   7.729 1.08e-14 ***
## Cigarette1Yes               -2.69414    0.16257 -16.572  < 2e-16 ***
## Marijuana1Yes               -2.77123    0.15199 -18.233  < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes    2.87373    0.16730  17.178  < 2e-16 ***
## Cigarette1Yes:Marijuana1Yes  3.22431    0.16098  20.029  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.461  on 7  degrees of freedom
## Residual deviance:   92.018  on 2  degrees of freedom
## AIC: 153.06
## 
## Number of Fisher Scoring iterations: 6

Pengujian Ada Tidaknya Interaksi antara X dan Z (Homogenous model vs Conditional Association on Y)

• Hipotesis

\[ H_0 : \lambda_{ik}^{XZ} = 0 \] Tidak ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Ganja (z)

\[ H_1 : \lambda_{ik}^{XZ} \ne 0 \] Ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Ganja (z)

• Tingkat Signifikansi

\[ \alpha = 5\% \]

• Statistik Uji

\[ \Delta \text{Deviance} = \text{Deviance model conditional on } Y - \text{Deviance model homogenous} \]

\[ = 92.018 - 0.372 = 91.644 \]

\[ db = \text{df model conditional on } Y - \text{df model homogenous} \]

\[ = 2 - 1 = 1 \]

• Daerah Penolakan

Tolak \(H_0\) jika \(\Delta\text{Deviance} > \chi^2_{0.05, 1} = 3.841\)

• Keputusan

\[ \text{Karena } 91.644 > 3.841, \text{ maka tolak } H_0 \] • Kesimpulan

Dengan taraf nyata 5%, ada interaksi yang signifikan antara Penggunaan Alkohol (X) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{ik}^{XZ}\).

Model 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_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

# Model Conditional Association on Z (Marijuana)
model_conditional_Z <- glm(counts ~ Alcohol + Cigarette + Marijuana +
              Alcohol:Marijuana + Cigarette:Marijuana,
              family = poisson(link = "log"), data = data)

summary(model_conditional_Z)
## 
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Marijuana + 
##     Cigarette:Marijuana, family = poisson(link = "log"), data = data)
## 
## Coefficients:
##                             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                  5.19207    0.06088  85.285  < 2e-16 ***
## Alcohol1Yes                  1.12719    0.06412  17.579  < 2e-16 ***
## Cigarette1Yes               -0.23512    0.05551  -4.235 2.28e-05 ***
## Marijuana1Yes               -6.62092    0.47370 -13.977  < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes    4.12509    0.45294   9.107  < 2e-16 ***
## Cigarette1Yes:Marijuana1Yes  3.22431    0.16098  20.029  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2851.46  on 7  degrees of freedom
## Residual deviance:  187.75  on 2  degrees of freedom
## AIC: 248.8
## 
## Number of Fisher Scoring iterations: 5

Pengujian Ada Tidaknya Interaksi antara X dan Y (Homogenous model vs Conditional Association on Z)

• Hipotesis

\[ H_0 : \lambda_{ij}^{XY} = 0 \] Tidak ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Rokok (y)

\[ H_1 : \lambda_{ij}^{XY} \ne 0 \] Ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Rokok (y)

• Tingkat Signifikansi

\[ \alpha = 5\% \]

• Statistik Uji

\[ \Delta\text{Deviance} = \text{Deviance model conditional on } Z - \text{Deviance model homogenous} \]

\[ = 187.75 - 0.374 = 187.380 \]

\[ db = \text{df model conditional on } Z - \text{df model homogenous} \]

\[ = 2 - 1 = 1 \]

• Daerah Penolakan

Tolak \(H_0\) jika \(\Delta\text{Deviance} > \chi^2_{0.05, 1} = 3.841\)

• Keputusan

\[ \text{Karena } 187.380 > 3.841, \text{ maka tolak } H_0 \] • Kesimpulan

Dengan taraf nyata 5%, Ada interaksi yang signifikan antara Penggunaan Alkohol (X) dan Penggunaan Rokok (Y)). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{ij}^{XY}\).

Pemilihan Model Terbaik

Ringkasan Model Log Linier

Model log linier dengan berbagai kombinasi parameter:

library(knitr)

model_table <- data.frame(
  Model = c(
    "Saturated",
    "Homogenous",
    "Conditional on X",
    "Conditional on Y",
    "Conditional on Z"
  ),
  Parameter = c(
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ", 
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ",             
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ",                        
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ",                        
    "λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ"                        
  ),
  Deviance = c(0.000, 0.374, 497.37, 92.018, 187.75),
  Param = c(8, 7, 6, 6, 6),
  df = c(0, 1, 1, 1, 1),
  AIC = c(65.043, 63.417, 558.41, 153.06, 248.80)
)

kable(
  model_table,
  format = "pipe",
  caption = "Ringkasan Model Log-Linear"
)
Ringkasan Model Log-Linear
Model Parameter Deviance Param df AIC
Saturated λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ 0.000 8 0 65.043
Homogenous λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ 0.374 7 1 63.417
Conditional on X λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ 497.370 6 1 558.410
Conditional on Y λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ 92.018 6 1 153.060
Conditional on Z λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ 187.750 6 1 248.800

Ringkasan Pengujian Interaksi 3 Arah dan 2 Arah

library(knitr)
uji_interaksi <- data.frame(
  Interaksi = c("XYZ", "YZ", "XZ", "XY"),
  Pengujian = c(
    "Saturated vs Homogenous",
    "Conditional on X vs Homogenous",
    "Conditional on Y vs Homogenous",
    "Conditional on Z vs Homogenous"
  ),
  "Δ Deviance" = c(0.374, 496.995, 91.644, 187.380),
  df = c(2, 1, 1, 1),
  ChiSqTabel = c(5.991, 3.841, 3.841, 3.841),
  Ket = c("Tidak ada interaksi", "Ada interaksi", "Ada interaksi", "Ada interaksi"),
  Keputusan = c("Terima $H_0$", "Tolak $H_0$", "Tolak $H_0$", "Tolak $H_0$")
)

# Tampilkan dengan kable
kable(uji_interaksi, format = "pipe", caption = "Ringkasan Uji Interaksi (Chi-Square)")
Ringkasan Uji Interaksi (Chi-Square)
Interaksi Pengujian Δ.Deviance df ChiSqTabel Ket Keputusan
XYZ Saturated vs Homogenous 0.374 2 5.991 Tidak ada interaksi Terima \(H_0\)
YZ Conditional on X vs Homogenous 496.995 1 3.841 Ada interaksi Tolak \(H_0\)
XZ Conditional on Y vs Homogenous 91.644 1 3.841 Ada interaksi Tolak \(H_0\)
XY Conditional on Z vs Homogenous 187.380 1 3.841 Ada interaksi Tolak \(H_0\)

Kesimpulan Pemilihan Model Terbaik

Dari hasil di atas diketahui bahwa asosiasi yang nyata terdapat pada model Homogenous dan ketiga model Conditional. Karena itulah, model akan dipilih berdasarkan nilai AIC yang paling kecil, dimana AIC terkecil terdapat pada model Homogenous yaitu sebesar 63.417. Oleh karena itu, model terbaik yaitu:

\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

Model terbaik adalah model log-linear tanpa interaksi tiga arah dan memuat semua interaksi dua arah.

Nilai Dugaan Model Terbaik

data.frame(
  Alcohol   = data$Alcohol,
  Cigarette = data$Cigarette,
  Marijuana = data$Marijuana,
  Frekuensi = data$Frekuensi,
  fitted    = fitted(model_homogenous)
)

Interpretasi:

Nilai fitted values merupakan prediksi dari model yang sudah mempertimbangkan interaksi dua arah. Nilai ini diperoleh dari estimasi model log-linear yang menyertakan semua efek utama dan semua interaksi dua arah. Dari perbandingan antara prediksi model dan nilai asli, didapatkan bahwa nilai tersebut saling mirip serta tidak terdapat perbedaan yang signifikan. Oleh karena itu, dapat disimpulkan bahwa model dapat memprediksi data dengan sangat baik.

Referensi

  • Agresti, A. (2002). Categorical data analysis (2nd ed.). Wiley-Interscience.

  • Agresti, A. (2013). Categorical data analysis (3rd ed.). John Wiley & Sons.

  • Agresti, A. (2019). An introduction to categorical data analysis (3rd ed.). Wiley.

  • Christensen, R. (1997). Log-linear models and logistic regression. Springer.

  • Fisher, R. A. (1935). The design of experiments. Oliver and Boyd.

  • Fox, J. (2008). Applied regression analysis and generalized linear models (2nd ed.). Sage.

  • Hilbe, J. M. (2009). Logistic regression models. Chapman and Hall/CRC.

  • Hosmer, D. W., Lemeshow, S., & Sturdivant, R. X. (2013). Applied logistic regression (3rd ed.). Wiley.

  • Howell, D. C. (2012). Statistical methods for psychology (8th ed.). Cengage Learning.

  • Kleinbaum, D. G., & Klein, M. (2010). Logistic regression: A self-learning text (3rd ed.). Springer.

  • Mehta, C. R., & Patel, N. R. (1983). A network algorithm for performing Fisher’s exact test in r × c contingency tables. Journal of the American Statistical Association, 78(382), 427–434.

  • Modul ADK, Pertemuan 13. (2025). Universitas Padjadjaran.

  • Powers, D. M. W. (2011). Evaluation: From precision, recall and F-measure to ROC, informedness, markedness and correlation. Journal of Machine Learning Technologies, 2(1), 37–63