Kata Pengantar

Puji dan syukur penulis panjatkan ke hadirat Tuhan Yang Maha Esa atas segala limpahan rahmat, karunia, dan petunjuk-Nya sehingga eBook ini yang berjudul Analisis Data Kategori dapat disusun dan diselesaikan dengan baik.

eBook ini disusun sebagai bentuk kontribusi dalam mendukung pembelajaran statistika terapan, khususnya dalam topik analisis data kategori. Penulis menyadari bahwa bidang ini memiliki peran penting dalam berbagai disiplin ilmu, terutama dalam memahami dan menganalisis data kualitatif yang sering dijumpai dalam praktik riset sosial, kesehatan, ekonomi, dan bidang bidang lainnya.

Penulis menyampaikan rasa terima kasih yang sebesar-besarnya kepada dosen mata kuliah Analisis data Kategori , Dr. I Gede Nyoman Mindra Jaya, atas segala ilmu, bimbingan, serta motivasi yang telah diberikan selama proses perkuliahan maupun dalam penyusunan materi eBook ini. Ulasan dan penjelasan beliau yang sangat membantu penulis dalam memahami dan menyusun konten secara lebih terarah.

Tujuan utama penulisan eBook ini adalah untuk menjadi referensi pembelajaran mandiri yang ringkas namun komprehensif, baik bagi mahasiswa, peneliti, maupun awam yang ingin memahami dasar-dasar serta aplikasi analisis data kategori. Penyusunan materi dilakukan dengan berusaha mengedepankan kejelasan konsep, contoh penerapan, serta relevansi dengan kebutuhan analisis nyata.

Akhir kata, penulis menyadari bahwa eBook ini masih jauh dari sempurna. Oleh karena itu, segala kritik dan saran yang membangun sangat diharapkan demi penyempurnaan karya ini di masa mendatang.

Semoga eBook ini dapat memberikan manfaat dan menjadi salah satu sumber belajar yang bermanfaat bagi siapa saja yang membutuhkannya.

1 Pendahuluan

1.1 Analisis Data Kategorik

Ilmu statistika telah berkembang sebagaimana perubahan dan perkembangan terjadi di berbagai disiplin ilmu. Dalam statistika, ukuran dan jenis - jenis data sangat beragam, salah satunya data kategorik. Data kategorik adalah jenis data yang mengelompokkan informasi ke dalam kelas-kelas diskrit, bukan angka kontinu.

Data ini mencakup berbagai label yang menjelaskan fenomena, seperti orientasi politik (progresif atau konservatif), selera hiburan (film, drama, atau reality show), dan segmentasi pengunjung website (baru atau kembali). Kategori-kategori ini tidak selalu memiliki urutan matematis, tapi kadang bisa jadi urut jika bersifat ordinal. Menganalisis data kategori penting untuk memahami pola, tren, dan hubungan antar variabel yang tidak hanya bisa dilihat dari angka. Dengan perkembangan machine learning dan big data, metode analisis yang melibatkan data kategori, seperti tabel kontingensi, uji chi-square, dan regresi logistik, terus diperbaiki untuk mendapatkan wawasan yang lebih dalam dan mendukung keputusan strategis.

1.2 Tujuan Utama Analisis Data Kategori

Analisis data kategori memiliki berbagai tujuan strategis yang mendukung riset dan praktik di sejumlah bidang. Beberapa tujuan utamanya meliputi:

  • Mengungkap Pola Tersembunyi

    Dengan mengelompokkan data ke dalam kategori, peneliti dapat menemukan kecenderungan yang tidak langsung terlihat, misalnya bagaimana preferensi digital berperilaku berbeda berdasarkan usia atau jenis kelamin.

  • Memahami Relasi antarvariabel

    Melalui uji statistik seperti chi-square, dapat diidentifikasi apakah terdapat hubungan yang signifikan antara dua atau lebih variabel, misalnya hubungan antara jenis langganan aplikasi musik digital dengan preferensi genre lagu.

  • Mendukung Pengambilan Keputusan

    Data kategori membantu pengambil keputusan untuk mengembangkan strategi yang lebih tepat sasaran. Contohnya, pemerintah dapat menggunakan informasi ini untuk menyusun program pelatihan kerja berdasarkan segmentasi pendidikan dan latar belakang sosial.

  • Membangun Model Prediktif

    Banyak algoritma prediktif, terutama regresi logistik, memanfaatkan data kategori untuk meramalkan kemungkinan kejadian, seperti memprediksi loyalitas pelanggan di sektor e-commerce berdasarkan data interaksi mereka.

1.3 Ruang Lingkup Data Kategori

Analisis data non-numerik ini mencakup dua tipe utama berdasarkan karakteristiknya, yaitu:

  • Nominal : Data yang tidak memiliki tingkat atau urutan, seperti jenis kendaraan (mobil, sepeda motor, sepeda) atau jenis kepribadian dalam studi psikologi.

  • Ordinal : Data yang memiliki urutan atau peringkat, contoh: tingkat kepuasan (sangat tidak puas, tidak puas, moderat, puas, sangat puas) atau level pendidikan (SD, SMP, SMA, Perguruan Tinggi).

Selain itu, data kategori dapat dibagi menjadi:

  • Biner : Hanya memiliki dua kategori, seperti sukses/gagal atau hadir/tidak hadir.

  • Multikategori: Melibatkan lebih dari dua kelompok, seperti klasifikasi jenis kejahatan (ringan, sedang, berat).

1.4 Aplikasi Data Kategori di Berbagai Bidang

Analisis data kategori tidak hanya penting di bidang akademik; penerapannya sudah merambah ke berbagai sektor yang semakin kompleks, misalnya:

1.4.1 Teknologi Digital dan Media Sosial

  • Contoh Kasus: Mengelompokkan pengguna berdasarkan jenis konten yang mereka sukai (video, artikel, podcast) dan menganalisis apakah preferensi tersebut berhubungan dengan tingkat interaksi di platform.

  • Aplikasi: Meningkatkan algoritma rekomendasi dengan menggunakan analisis segmentasi berbasis kategori perilaku pengguna.

1.4.2 Kesehatan dan Kesejahteraan

  • Contoh Kasus: Penelitian tentang pengaruh pola makan (vegetarian, vegan, omnivora) terhadap kejadian penyakit kardiovaskular.

  • Aplikasi: Membantu klinik menentukan program pencegahan yang sesuai berdasarkan kategori gaya hidup pasien.

1.4.3 Pemasaran dan E-Commerce

  • Contoh Kasus: Segmentasi pelanggan berdasarkan perilaku pembelian (impulsif/terencana) dan preferensi produk (fashion, elektronik, makanan).

  • Aplikasi: Pengolahan data tersebut memungkinkan perusahaan menyesuaikan strategi promosi dan meningkatkan inferensi model prediktif mengenai loyalitas pelanggan.

1.4.4 Analisis Perilaku Konsumen di Era Digital

  • Contoh Kasus: Survei interaktif yang mengelompokkan responden berdasarkan sikap terhadap inovasi teknologi (antusias, skeptis, netral) dan menghubungkannya dengan adopsi produk terbaru.

  • Aplikasi: Membantu memprediksi tren pasar dan memberikan insight pada pengembangan produk baru.

1.4.5 Kebijakan Publik dan Pemerintahan

  • Contoh Kasus: Analisis tanggapan masyarakat terhadap kebijakan ramah lingkungan, dikategorikan sebagai mendukung, ragu, atau menolak.

  • Aplikasi: Menjadi dasar bagi perancang kebijakan dalam menyesuaikan inisiatif publik sehingga lebih responsif terhadap aspirasi warga.

2 Metode Analisis Data Kategori

Bergantung pada tujuan riset, terdapat berbagai metode analisis yang spesifik untuk data kategori, antara lain:

  • Tabel Kontingensi dan Uji Chi-Square: Teknik dasar yang digunakan untuk menguji apakah distribusi kategori antar variabel bersifat independen atau memiliki hubungan, misalnya analisis hubungan antara jenis konten dan interaksi pengguna.

  • Regresi Logistik: Digunakan untuk memodelkan probabilitas kejadian biner, contohnya dalam meramalkan kemungkinan terjadinya pembelian ulang pelanggan.

  • Analisis Correspondence: Teknik ini memvisualisasikan hubungan kompleks antar kategori, seperti asosiasi antara jenis hiburan dan preferensi demografis.

  • Metode Pembelajaran Mesin: Algoritma seperti Decision Tree dan Random Forest mampu mengklasifikasikan data kategori dengan tingkat akurasi tinggi, contohnya dalam mengidentifikasi segmen pasar potensial berdasarkan aktivitas digital.

3 Distribusi Peluang dalam Data Kategori

Distribusi peluang dalam analisis data kategori merupakan dasar untuk memahami bagaimana kemungkinan suatu kategori muncul dalam suatu sampel atau populasi. Beberapa distribusi probabilitas yang sering digunakan meliputi Bernoulli, Binomial, Multinomial, dan Poisson, yang masing-masing memiliki peran penting dalam pemodelan kejadian diskrit.

3.1 Distribusi Bernoulli

Distribusi Bernoulli digunakan untuk kejadian biner (ya/tidak, sukses/gagal, benar/salah). Jika suatu peristiwa memiliki probabilitas sukses sebesar 𝑝 , maka probabilitas kegagalan adalah 1− 𝑝 .

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

Dengan :

\(p\) adalah probabilitas sukses (\(0 \leq p \leq 1\))

\(1-p\) adalah probabilitas gagal

\(k\) adalah banyaknya kejadian sukses

Simulasi R:

# Simulasi Distribusi Bernoulli (10 percobaan, p=0.5)
set.seed(123)
bernoulli_sample <- rbinom(n = 10, size = 1, prob = 0.5)
bernoulli_sample
##  [1] 0 1 0 1 1 0 1 1 1 0

3.2 Distribusi Binomial

Distribusi Binomial adalah perluasan dari Bernoulli untuk kasus lebih dari satu percobaan independen. Jika dalam 𝑛 percobaan, setiap percobaan memiliki probabilitas sukses 𝑝 , maka jumlah sukses dalam percobaan mengikuti distribusi binomial.

Rumus:

\(P(X = k)\) = \(\binom{n}{k} p^k (1-p)^{n-k}\) Dengan :

  • \(n\) adalah jumlah percobaan,

  • \(k\) adalah jumlah keberhasilan,

  • \(p\) adalah probabilitas keberhasilan dalam satu percobaan,

  • \((1-p)\) adalah probabilitas kegagalan.

Simulasi R:

# Simulasi Distribusi Binomial (25 percobaan, n = 10, probabilitas sukses p = 0.5)
set.seed(123)
binomial_sample <- rbinom(n = 25, size = 10, prob = 0.5)
binomial_sample
##  [1] 4 6 5 7 7 2 5 7 5 5 8 5 6 5 3 7 4 2 4 8 7 6 6 9 6

3.3 Distribusi Multinomial

Distribusi Multinomial generalisasi dari distribusi binomial untuk kasus dengan lebih dari dua hasil (multiple outcomes). Ini digunakan ketika suatu eksperimen memiliki k ≥ 2 kemungkinan hasil, masing-masing dengan probabilitas tetap.

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

\(n\) : Banyaknya percobaan.

\(k\) : Jumlah kategori/kelas yang mungkin.

\(xᵢ\) : Banyaknya kejadian kategori ke-i.

\(pᵢ\): Probabilitas kategori ke-i (harus memenuhi ∑pᵢ = 1).

Simulasi R:

# Simulasi Distribusi Multinomial (1 percobaan, 3 kategori, total percobaan size = 15)
set.seed(123)
# rmultinom menghasilkan matriks dengan masing-masing kolom mewakili satu simulasi
multinomial_sample <- rmultinom(n = 1, size = 15, prob = c(0.3, 0.5, 0.2))
multinomial_sample
##      [,1]
## [1,]    3
## [2,]    7
## [3,]    5

3.4 Distribusi Poisson

Distribusi Poisson sering digunakan untuk memodelkan jumlah kejadian dalam interval waktu atau ruang tertentu.

Rumus: \[P(X = k) = \frac{\lambda^k e^{-\lambda}}{k!}\]

Dengan :

  • \(k\) adalah jumlah kejadian (non-negative integer),

  • \(\lambda\) adalah Rata-rata kejadian dalam interval

  • \(e\) adalah bilangan Euler’s (~2.71828).

Simulasi R:

# Parameter
lambda <- 4    # average rate (λ)
k <- 8         # number of events

# P(X = k)
sampel <- rpois(k, lambda)
print(sampel)
## [1] 3 6 7 1 4 7 4 4

4 Desain Sampling Data Kategori

Dalam penelitian kategori dan epidemiologi, pemilihan desain sampling sangat krusial untuk menghasilkan inferensi yang valid. Dua pendekatan utama adalah prospective sampling dan retrospective sampling. Setiap pendekatan memiliki kekuatan dan keterbatasan tersendiri, serta metode pengambilan sampel yang berbeda sesuai dengan tujuan studi.

4. 1 Prospective Sampling

Prospective sampling (atau desain prospektif) adalah metode di mana subjek penelitian dipilih sebelum outcome (hasil) terjadi. Data dikumpulkan secara luring (forward-looking), sehingga memungkinkan peneliti mengikuti perkembangan subjek dari waktu pemilihan sampel ke masa depan.

Jenis-jenis Prospective Sampling

Beberapa jenis desain prospective sampling antara lain:

  • Cohort Study (Studi Kohort):
    Subjek dipilih berdasarkan paparan tertentu dan diikuti secara longitudinal untuk mengamati kejadian outcome.
    Contoh: Mengikuti sekelompok perokok dan non-perokok selama 10 tahun untuk melihat kejadian kanker paru-paru.

  • Randomized Controlled Trial (RCT):
    Subjek diacak ke dalam grup intervensi dan kontrol, kemudian outcome diukur ke depannya.
    Contoh: Menguji efektivitas vaksin dengan mengacak peserta ke grup vaksin dan placebo.

  • Prospective Convenience Sampling:
    Sampling dilakukan secara prospektif berdasarkan ketersediaan subjek—meskipun kurang representatif, metode ini cepat dan mudah dijalankan dalam keadaan terbatas.

4.2 Retrospective Sampling

Retrospective sampling (atau desain retrospektif) merupakan metode di mana data dikumpulkan setelah outcome telah terjadi. Pendekatan ini umumnya menggunakan data yang sudah ada (existing records) untuk mengidentifikasi subjek berdasarkan status outcome lalu menelusuri faktor paparan secara mundur.

Jenis-jenis Retrospective Sampling

Beberapa jenis desain retrospective sampling antara lain:

  • Case-Control Study (Studi Kasus-Kontrol):
    Subjek dipilih berdasarkan outcome (misalnya, penyakit) dan kontrol tanpa outcome, kemudian informasi mengenai paparan dikumpulkan secara retrospektif.
    Contoh: Mengidentifikasi pasien dengan kanker paru-paru (kasus) dan tanpa kanker paru-paru (kontrol) serta membandingkan paparan terhadap rokok di masa lampau.

  • Retrospective Cohort Study:
    Data historis digunakan untuk mengelompokkan subjek berdasarkan paparan yang terjadi di masa lalu dan kemudian mengamati outcome yang terjadi setelahnya.
    Contoh: Menganalisis arsip rumah sakit untuk menentukan tingkat kematian pada pasien dengan paparan tertentu.

  • Cross-Sectional Retrospective Study:
    Pengumpulan data dilakukan pada satu titik waktu, tetapi dengan informasi tentang sejarah paparan dan outcome yang dikumpulkan sekaligus.

5 Tabel Kontingensi 2x2

Dalam aplikasi analisis data kategori, yang dimaksud dengan tabel kontingensi adalah tabel yang menampilkan distribusi frekuensi dari dua atau lebih variabel kategorik secara bersamaan. Tabel ini mengorganisasikan data ke dalam baris dan kolom, di mana masing-masing baris dan kolom mewakili kategori dari variabel yang berbeda. Dengan demikian, tabel kontingensi memudahkan analisis hubungan atau asosiasi antara variabel-variabel tersebut.

Tabel kontingensi 2x2 akan mengorganisasikan data dalam empat sel:

Struktur Tabel Kontingensi 2x2
Kategori 1 (+) Kategori 2 (-) Jumlah
Kelompok 1 a b a + b
Kelompok 2 c d c + d
Jumlah a + c b + d N

Contoh kasus kontingensi 2x2

Hubungan Pemberian ASI dengan Kejadian Gizi Kurang

Contoh kasus kontingensi 2x2

Kondisi Gizi

Pemberian ASI

Gizi Kurang Gizi Baik Jumlah
ASI Ekslusif 9 23 32
ASI Non Ekslusif 21 7 28
Jumlah 30 30 60

Pada kasus di atas, terdapat data yang menggambarkan ter-jadinya kejadian kekurangan gizi atau umumnya gizi buruk berdasarkan pemberian ASI (Air Susu Ibu). Terdapat dua kelompok pada penelitian ini yaitu kelompok ASI Ekslusif () dan ASI Non Ekslusif, adapun arti ekslusif dalam konteks kasus ini adalah pemberian ASI …. . Sedangkan kejadian yang diamati yaitu Gizi Kurang dan Gizi Baik.

5.1 Distribusi Peluang Tabel Kontingensi 2x2

5.1.1 Peluang Bersama

Peluang bersama dapat diartikan sebagai peluang terjadinya dua kejadian pada satu sel tabel secara bersamaan. Peluang bersama dihitung dengan rumus :

\[ P(A\cap B) = \frac{n_{ij}}{N} \]

Contoh Perhitungan :

Peluang bersama dihitung dengan membagi frekuensi tiap sel dengan total sampel (N = 60).

  • P(Ekslusif, Gizi Kurang): \[ P(\text{Ekslusif, Gizi Kurang}) = \frac{9}{60} = 0.15 \]

  • P(Ekslusif, Gizi Baik): \[ P(\text{Ekslusif, Gizi Baik}) = \frac{23}{60} \approx 0.3833 \]

  • P(Non Eksklusif, Gizi Kurang): \[ P(\text{Non Ekslusif, Gizi Kurang}) = \frac{21}{60} = 0.35 \]

  • P(Non Eksklusif, Gizi Baik): \[ P(\text{Non Ekslusif, Gizi Baik}) = \frac{7}{60} \approx 0.1167 \]

# Membuat matriks data sesuai tabel
data_observed <- matrix(c(9, 23, 21, 7), nrow = 2, byrow = TRUE)
rownames(data_observed) <- c("ASI Ekslusif", "ASI Non Ekslusif")
colnames(data_observed) <- c("Gizi Kurang", "Gizi Baik")
data_observed
##                  Gizi Kurang Gizi Baik
## ASI Ekslusif               9        23
## ASI Non Ekslusif          21         7
total <- sum(data_observed) 

# Menghitung peluang bersama tiap sel
joint_prob <- data_observed / total
joint_prob
##                  Gizi Kurang Gizi Baik
## ASI Ekslusif            0.15 0.3833333
## ASI Non Ekslusif        0.35 0.1166667
total <- sum(data_observed)  # Total sampel = 60
joint_prob <- data_observed / total
print(joint_prob)
##                  Gizi Kurang Gizi Baik
## ASI Ekslusif            0.15 0.3833333
## ASI Non Ekslusif        0.35 0.1166667
# Ekstraksi nilai menggunakan nama baris dan kolom yang terverifikasi
P_EKS_GK <- joint_prob["ASI Ekslusif", "Gizi Kurang"]     # 9/60
P_EKS_GB <- joint_prob["ASI Ekslusif", "Gizi Baik"]        # 23/60
P_NON_GK <- joint_prob["ASI Non Ekslusif", "Gizi Kurang"]   # 21/60
P_NON_GB <- joint_prob["ASI Non Ekslusif", "Gizi Baik"]     # 7/60

cat("P(Eksklusif, Gizi Kurang) =", P_EKS_GK, "\n")
## P(Eksklusif, Gizi Kurang) = 0.15
cat("P(Eksklusif, Gizi Baik) =", P_EKS_GB, "\n")
## P(Eksklusif, Gizi Baik) = 0.3833333
cat("P(Non Ekslusif, Gizi Kurang) =", P_NON_GK, "\n")
## P(Non Ekslusif, Gizi Kurang) = 0.35
cat("P(Non Ekslusif, Gizi Baik) =", P_NON_GB, "\n")
## P(Non Ekslusif, Gizi Baik) = 0.1166667

5.1.2 Peluang Marginal

Peluang terjadinya sebuah kejadian (event) tanpa mempertimbangkan kejadian lainnya. Peluang ini diperoleh dengan menjumlahkan probabilitas gabungan (joint probability) pada seluruh nilai variabel lainnya. Dihitung dengan rumus :

\[ P(A) = \frac{a+b}{N}, \quad P(B) = \frac{a+c}{N} \]

Contoh perhitungan :

  • P(Ekslusif):
    \[ P(\text{Ekslusif}) = \frac{32}{60} \approx 0.5333 \]

  • P(Non Ekslusif):
    \[ P(\text{Non Ekslusif}) = \frac{28}{60} \approx 0.4667 \]

  • P(Gizi Kurang):
    \[ P(\text{Gizi Kurang}) = \frac{30}{60} = 0.5 \]

  • P(Gizi Baik):
    \[ P(\text{Gizi Baik}) = \frac{30}{60} = 0.5 \]

Dalam R :

# Total per baris (jenis ASI) dan per kolom (status gizi)
marginal_rows <- margin.table(data_observed, 1)  
marginal_cols <- margin.table(data_observed, 2)  

# Peluang marginal
p_exklusif <- marginal_rows["ASI Ekslusif"] / total
p_non_eksklusif <- marginal_rows["ASI Non Ekslusif"] / total
p_gizi_kurang  <- marginal_cols["Gizi Kurang"] / total
p_gizi_baik    <- marginal_cols["Gizi Baik"] / total

cat("P(Eksklusif) =", p_exklusif, "\n")
## P(Eksklusif) = 0.5333333
cat("P(Non Eksklusif) =", p_non_eksklusif, "\n")
## P(Non Eksklusif) = 0.4666667
cat("P(Gizi Buruk / Gizi Kurang) =", p_gizi_kurang, "\n")
## P(Gizi Buruk / Gizi Kurang) = 0.5
cat("P(Gizi Baik) =", p_gizi_baik, "\n")
## P(Gizi Baik) = 0.5

5.1.3 Peluang Bersyarat

Peluang bersyarat menghitung peluang terjadinya suatu kejadian yang terikat pada kejadian lain yang sudah terjadi. Dihitung dengan rumus :

\[ P(B | A) = \frac{P(A \cap B)}{P(A)} \]

Contoh perhitungannya sebagai berikut :

  • P(Gizi Kurang | Non Ekslusif):
    \[ P(\text{Gizi Kurang} \mid \text{Non Ekslusif}) = \frac{21}{28} = 0.75 \]

  • P(Gizi Kurang | Ekslusif):
    \[ P(\text{Gizi Kurang} \mid \text{Ekslusif}) = \frac{9}{32} \approx 0.2813 \]

Dalam R :

# Peluang Gizi Kurang berdasarkan jenis ASI
cond_non_ekslusif <- data_observed["ASI Non Ekslusif", "Gizi Kurang"] /
                      marginal_rows["ASI Non Ekslusif"]
cond_ekslusif <- data_observed["ASI Ekslusif", "Gizi Kurang"] /
                 marginal_rows["ASI Ekslusif"]

cat("P(Gizi Kurang | Non Ekslusif) =", cond_non_ekslusif, "\n")
## P(Gizi Kurang | Non Ekslusif) = 0.75
cat("P(Gizi Kurang | Eksklusif) =", cond_ekslusif, "\n")
## P(Gizi Kurang | Eksklusif) = 0.28125

5.2 Ukuran Asosiasi Tabel Kontingensi 2x2

Ukuran asosiasi bertujuan untuk menentukan hubungan antar variabel, dalam kasus kontingensi 2x2 asosiasi menggambarkan hubungan statistik antara dua variabel dan dua kejadian.

Adapun ukuran asosiasi yang umum digunakan pada kasus kontingensi 2x2 yaitu :

5.2.1 Risk Difference (RD)

Mengukur selisih risiko antar kelompok atau menggambarkan perbedaan peluang kejadian suatu hasil pada kelompok yang berbeda

\[ RD = P(B|A) - P(B|\neg A) \]

  • Jika 𝑅𝐷>0, maka risiko kejadian lebih tinggi di Grup 1 dibandingkan Grup 2.

  • Jika 𝑅𝐷<0, maka risiko kejadian lebih rendah di Grup 1 dibandingkan Grup 2.

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

Untuk ASI Ekslusif:

\[ R_{\text{Ekslusif}} = \frac{9}{32} \approx 0.2813. \]

Untuk ASI Non Ekslusif:

\[ R_{\text{Non Ekslusif}} = \frac{21}{28} = 0.75. \]

Perhitungan R:

# Data:
a1 <- 9     
n1 <- 32   
a2 <- 21    
n2 <- 28    

# Risiko
risk_exklusif <- a1 / n1
risk_non_ekslusif <- a2 / n2

cat("Risk (ASI Ekslusif) =", risk_exklusif, "\n")
## Risk (ASI Ekslusif) = 0.28125
cat("Risk (ASI Non Ekslusif) =", risk_non_ekslusif, "\n")
## Risk (ASI Non Ekslusif) = 0.75
# Risk Difference (RD)
RD <- risk_non_ekslusif - risk_exklusif
cat("Risk Difference (RD) =", RD, "\n")
## Risk Difference (RD) = 0.46875

Dari hasil perhitungan RD, dapat dilihat bahwa RD Non Ekslusif > RD Eklsusif yang artinya risiko kejadian gizi kurang pada kelompok ASI Non Ekslusif lebih tinggi dibandingkan dengan risiko kejadian gizi kurang pada kelompok ASI Eklsusif.

5.2.2 Relative Risk (RR)

\[ RR = \frac{(n11 / n1.)}{(n21/n2.)} \]

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

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

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

\[ RR = \frac{R_{\text{Non Ekslusif}}}{R_{\text{Ekslusif}}} \approx \frac{0.75}{0.2813} \approx 2.6667. \]

Perhitungan R:

# Relative Risk (RR)
RR <- risk_non_ekslusif / risk_exklusif
cat("Relative Risk (RR) =", RR, "\n")
## Relative Risk (RR) = 2.666667

Berdasarkan nilai RR yang didapatkan, disimpulkan bahwa frekuensi terjadinya kejadian gizi kurang lebih sering terjadi pada kelompok ASI Non ekslusif.

5.2.3 Odds Ratio (OR)

\[ OR = \frac{(a/c)}{(b/d)} = \frac{a \times d}{b \times c} \]

Perhitungan Odds pada masing - masing kelompok:

  • Odds pada ASI Ekslusif:

\[ \text{Odds}_{\text{Ekslusif}} = \frac{9}{23} \approx 0.3913. \]

  • Odds pada ASI Non Ekslusif:

\[ \text{Odds}_{\text{Non Ekslusif}} = \frac{21}{7} = 3. \]

Kemudian,

\[ OR = \frac{\text{Odds}_{\text{Non Ekslusif}}}{\text{Odds}_{\text{Ekslusif}}} \approx \frac{3}{0.3913} \approx 7.6667. \]

Perhitungan dalam R:

# Menghitung odds di setiap kelompok
odds_exklusif <- a1 / (n1 - a1)  # odds = 9 / (32 - 9) = 9/23
odds_non_ekslusif <- a2 / (n2 - a2)  # odds = 21 / (28 - 21) = 21/7 = 3

cat("Odds (ASI Ekslusif) =", odds_exklusif, "\n")
## Odds (ASI Ekslusif) = 0.3913043
cat("Odds (ASI Non Ekslusif) =", odds_non_ekslusif, "\n")
## Odds (ASI Non Ekslusif) = 3
# Odds Ratio (OR)
OR <- odds_non_ekslusif / odds_exklusif
cat("Odds Ratio (OR) =", OR, "\n")
## Odds Ratio (OR) = 7.666667

Nilai OR sebesar 7,6667 menunjukkan bahwa peluang kejadian gizi kurang pada kelompok ASI Non eklsusif lebih besar dibandingkan kelompok ASI Eklsusif.

6 Inferensi Tabel Kontingensi 2 Arah

Pengambilan kesimpulan berdasarkan hasil analisis sampel dari populasi umumnya disebut sebagai inferensi. Pengaplikasian inferensi pada kasus tabel kontingensi 2 arah mengacu pada analisis hubungan dua variabel kategorik pada tabel kontingensi. Inferensi ini bertujuan untuk memahami hubungan antara variabel - variabel yang bersangkutan dengan cara menghitung estimasi dan menguji hipotesis.

6.1 Estimasi

Estimasi adalah proses menaksir parameter populasi menggunakan data sampel. Estimasi terbagi menjadi dua jenis:

  • Estimasi Titik: Memberikan satu nilai sebagai perkiraan terbaik untuk parameter.
  • Estimasi Interval: Menyediakan rentang nilai yang diyakini mencakup parameter populasi sebenarnya, dengan tingkat kepercayaan tertentu, misalnya 95%.

Di sini, kita akan membahas dua pendekatan estimasi utama: estimasi titik (point estimation) dan estimasi interval (interval estimation). Untuk estimasi titik pada tabel kontingensi dua arah, proporsi biasanya merujuk pada proporsi sel tertentu atau kombinasi sel.

6.1.1 Estimasi Titik

Estimasi titik digunakan untuk menentukan satu nilai spesifik sebagai perkiraan terbaik dari parameter populasi. Sebagai contoh, untuk memperkirakan proporsi populasi dari data kategorik, kita dapat menggunakan:

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

di mana:

\(x\) adalah jumlah pengamatan kejadian, dan \(n\) adalah ukuran sampel.

6.1.2 Estimasi Interval

Estimasi interval memberikan rentang nilai yang diyakini mengandung parameter populasi dengan tingkat keyakinan tertentu (contoh: 95%).
Untuk proporsi, interval kepercayaan dapat dihitung dengan:

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

di mana: - \(\hat{p}\) adalah estimasi proporsi, - \(z_{\alpha/2}\) adalah nilai kritis dari distribusi normal untuk tingkat kepercayaan yang diinginkan, - \(n\) adalah ukuran sampel.

6.2 Uji Hipotesis

6.2.1 Uji Proporsi

Uji proporsi digunakan untuk mengevaluasi perbedaan proporsi kejadian antara dua kelompok yang disajikan dalam tabel kontingensi. Metode ini bertujuan untuk menentukan apakah terdapat perbedaan signifikan dalam proporsi kejadian antara kedua kelompok tersebut.
Statistik uji untuk dua proporsi dapat dihitung dengan:

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

di mana: - \(\hat{p}_1\) dan \(\hat{p}_2\) adalah estimasi proporsi masing-masing sampel, - \(p\) adalah proporsi gabungan dari kedua sampel, - \(n_1\) dan \(n_2\) adalah ukuran masing-masing sampel.

Untuk menggambarkan dan menelusuri uji proporsi lebih jelas, akan dilakukan perhitungan dengan data berikut :

Kondisi Gizi

Pemberian ASI

Gizi Kurang Gizi Baik Jumlah
ASI Ekslusif 9 23 32
ASI Non Ekslusif 21 7 28
Jumlah 30 30 60

Hipotesis

H0 : Tidak terdapat perbedaan proporsi antar kelompok

H1 : Terdapat perbedaan proporsi antar kelompok

Estimasi Proporsi :

  • Proporsi pada ASI Ekslusif:

\[ p_1 = \frac{9}{32} \approx 0.2813. \]

  • Proporsi pada ASI Non Ekslusif:

\[ p_2 = \frac{21}{28} \approx 0.75. \]

Proporsi Gabungan:

Hitung proporsi gabungan:

\[ p = \frac{9 + 21}{32 + 28} = \frac{30}{60} = 0.5. \]

Standar Error (SE):

Hitung standar error menggunakan formula berikut:

\[ SE = \sqrt{p(1-p)\left(\frac{1}{n_1}+\frac{1}{n_2}\right)}. \]

Substitusi nilai:

\[ SE = \sqrt{0.5 \times 0.5 \times \left(\frac{1}{32} + \frac{1}{28}\right)}. \]

\[ SE \approx 0.1294. \]

Statistik Uji Z:

Hitung statistik uji Z:

\[ Z = \frac{p_1 - p_2}{SE}. \]

Substitusi nilai:

\[ Z = \frac{0.2813 - 0.75}{0.1294} \approx -3.62. \]

Perhitungan R :

set.seed(123)
data_matrix <- matrix(c(50, 30, 30, 50), nrow = 2, byrow = TRUE)
dimnames(data_matrix) <- list("Terpapar" = c("Ya", "Tidak"),
                              "Kejadian" = c("Ya", "Tidak"))
print(data_matrix)
##         Kejadian
## Terpapar Ya Tidak
##    Ya    50    30
##    Tidak 30    50
# Uji Proporsi dengan variabel yang benar
# x: jumlah kejadian "Ya" pada kolom "Kejadian" untuk masing-masing baris (Terpapar)
# n: total masing-masing baris dari data_matrix
prop_test <- prop.test(x = c(data_matrix["Ya", "Ya"], data_matrix["Tidak", "Ya"]),
                       n = c(sum(data_matrix["Ya", ]), sum(data_matrix["Tidak", ])))
print(prop_test)
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(data_matrix["Ya", "Ya"], data_matrix["Tidak", "Ya"]) out of c(sum(data_matrix["Ya", ]), sum(data_matrix["Tidak", ]))
## X-squared = 9.025, df = 1, p-value = 0.002663
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  0.08747151 0.41252849
## sample estimates:
## prop 1 prop 2 
##  0.625  0.375

Dari hasil perhitungan, didapatkan nilai Z sebesar -3,62 dan nilai p - value sebesar 0,0026 yang bilamana jika nilai p - value dibandingkan dengan taraf signifikansi 0,05 maka didapatkan hasil yang signifikan dan dapat diambil keputusan untuk tolak H0 , sehingga dapat disimpulkan bahwa proporsi kedua kelompok berbeda signifikan.

6.2.2 Uji Asosiasi

Uji asosiasi pada kasus kontingensi 2 arah berguna untuk mengukur hubungan antara dua variabel kategorik. Seperti yang sudah dijelaskan pada bagian sebelumnya, ukuran asosiasi terbagi menjadi 3 yaitu :

Sedangkan perumusan hipotesis untuk pengujian hipotesis pada ukuran asosiasi yaitu :

H0 : Tidak ada asosiasi antar dua variabel

H1 : Terdapat asosiasi antar dua variabel

Perhitungan Manual

Hitung Risk Difference:

\[ \text{Risk Difference} = \text{Risiko}_{\text{Eksklusif}} - \text{Risiko}_{\text{Non Eksklusif}} \]

\[ \text{Risiko}_{\text{Eksklusif}} = \frac{9}{32} \approx 0.2813 \]

\[ \text{Risiko}_{\text{Non Eksklusif}} = \frac{21}{28} \approx 0.75 \]

\[ \text{Risk Difference} = 0.2813 - 0.75 = -0.4687 \]

Hitung Standard Error Risk Difference:

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

Substitusi nilai:

\[ SE_{\text{RD}} = \sqrt{\frac{0.2813 \times 0.7187}{32} + \frac{0.75 \times 0.25}{28}} \approx \sqrt{0.00633 + 0.00669} \approx \sqrt{0.01302} \approx 0.1141 \]

Hitung Statistik Uji Z untuk Risk Difference:

\[ Z_{\text{RD}} = \frac{\text{Risk Difference}}{SE_{\text{RD}}} \]

\[ Z_{\text{RD}} = \frac{-0.4687}{0.1141} \approx -4.11 \]

1.5 Relative Risk

Hitung Relative Risk:

\[ \text{Relative Risk} = \frac{\text{Risiko}_{\text{Eksklusif}}}{\text{Risiko}_{\text{Non Eksklusif}}} \]

\[ \text{Relative Risk} = \frac{0.2813}{0.75} \approx 0.3751 \]

Hitung Standard Error untuk Log(Relative Risk): \[ SE_{\log(\text{RR})} = \sqrt{\frac{1-p_1}{p_1n_1} + \frac{1-p_2}{p_2n_2}} \]

Substitusi nilai: \[ SE_{\log(\text{RR})} = \sqrt{\frac{1-0.2813}{0.2813 \times 32} + \frac{1-0.75}{0.75 \times 28}} \approx \sqrt{0.0774 + 0.0476} \approx \sqrt{0.125} \approx 0.3536 \]

Hitung Statistik Uji Z untuk Relative Risk:

\[ Z_{\text{RR}} = \frac{\log(\text{Relative Risk})}{SE_{\log(\text{RR})}} \]

\[ \log(\text{Relative Risk}) = \log(0.3751) \approx -0.9803 \]

\[ Z_{\text{RR}} = \frac{-0.9803}{0.3536} \approx -2.77 \]

1.6 Odds Ratio (OR) Hitung Odds Ratio:

\[ \text{Odds Ratio} = \frac{\text{Odds}_{\text{Eksklusif}}}{\text{Odds}_{\text{Non Eksklusif}}} \]

\[ \text{Odds}_{\text{Eksklusif}} = \frac{9}{23} \approx 0.3913 \]

\[ \text{Odds}_{\text{Non Eksklusif}} = \frac{21}{7} = 3 \]

\[ \text{Odds Ratio} = \frac{0.3913}{3} \approx 0.1304 \]

Hitung Standard Error untuk Log(Odds Ratio):

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

Dengan \(a = 9\), \(b = 23\), \(c = 21\), \(d = 7\):

\[ SE_{\log(\text{OR})} = \sqrt{\frac{1}{9} + \frac{1}{23} + \frac{1}{21} + \frac{1}{7}} \approx \sqrt{0.1111 + 0.0435 + 0.0476 + 0.1429} \approx \sqrt{0.3451} \approx 0.5872 \]

Hitung Statistik Uji Z untuk Odds Ratio:

\[ Z_{\text{OR}} = \frac{\log(\text{Odds Ratio})}{SE_{\log(\text{OR})}} \]

\[ \log(\text{Odds Ratio}) = \log(0.1304) \approx -2.0380 \]

\[ Z_{\text{OR}} = \frac{-2.0380}{0.5872} \approx -3.47 \] Perhitungan R :

# Data awal
success1 <- 9
n1 <- 32
success2 <- 21
n2 <- 28
a <- success1
b <- n1 - success1
c <- success2
d <- n2 - success2

# Risk Difference
p1 <- success1 / n1
p2 <- success2 / n2
risk_difference <- p1 - p2
se_rd <- sqrt((p1 * (1 - p1)) / n1 + (p2 * (1 - p2)) / n2)
z_rd <- risk_difference / se_rd

# Relative Risk
relative_risk <- p1 / p2
se_log_rr <- sqrt((1 - p1) / (p1 * n1) + (1 - p2) / (p2 * n2))
log_rr <- log(relative_risk)
z_rr <- log_rr / se_log_rr

# Odds Ratio
odds1 <- a / b
odds2 <- c / d
odds_ratio <- odds1 / odds2
se_log_or <- sqrt(1/a + 1/b + 1/c + 1/d)
log_or <- log(odds_ratio)
z_or <- log_or / se_log_or

# Tampilkan hasil
cat("Risk Difference:", risk_difference, "\n")
## Risk Difference: -0.46875
cat("SE Risk Difference:", se_rd, "\n")
## SE Risk Difference: 0.114077
cat("Z Risk Difference:", z_rd, "\n\n")
## Z Risk Difference: -4.109066
cat("Relative Risk:", relative_risk, "\n")
## Relative Risk: 0.375
cat("SE Log(Relative Risk):", se_log_rr, "\n")
## SE Log(Relative Risk): 0.3029288
cat("Z Relative Risk:", z_rr, "\n\n")
## Z Relative Risk: -3.237821
cat("Odds Ratio:", odds_ratio, "\n")
## Odds Ratio: 0.1304348
cat("SE Log(Odds Ratio):", se_log_or, "\n")
## SE Log(Odds Ratio): 0.5874228
cat("Z Odds Ratio:", z_or, "\n")
## Z Odds Ratio: -3.467489

Kesimpulan:

  • Risk Difference (RD) mencerminkan selisih risiko absolut antara dua kelompok.

  • Relative Risk (RR) menunjukkan perbandingan peluang terjadinya suatu kejadian pada kedua kelompok.

  • Odds Ratio (OR) mengevaluasi perbandingan peluang suatu kejadian antara dua kelompok.

  • Standard error dan uji statistik Z digunakan untuk menguji apakah ukuran asosiasi ini signifikan secara statistik.

6.2.3 Uji Independensi

Uji Independensi bertujuan untuk menentukan apakah ada hubungan antara dua variabel kategorik dalam tabel kontingensi. Uji chi-square adalah salah satu uji yang paling umum, dengan statistik uji:

\[ \chi^2 = \sum_{i,j} \frac{(O_{ij} - E_{ij})^2}{E_{ij}} \]

di mana: - \(O_{ij}\) adalah frekuensi yang diobservasi pada sel \((i,j)\), - \(E_{ij}\) adalah frekuensi yang diharapkan dihitung dengan:

\[ E_{ij} = \frac{(\text{Total baris } i) \times (\text{Total kolom } j)}{\text{Total keseluruhan}} \] 1. Expected Counts (E)

Expected count untuk setiap sel dihitung dengan rumus:

\[ E_{ij} = \frac{(\text{Row Total}_i)(\text{Column Total}_j)}{\text{Grand Total}} \]

  • Untuk ASI Eksklusif dan Gizi Kurang: \[ E_{11} = \frac{32 \times 30}{60} = 16 \]
  • Untuk ASI Eksklusif dan Gizi Baik: \[ E_{12} = \frac{32 \times 30}{60} = 16 \]
  • Untuk ASI Non Eksklusif dan Gizi Kurang: \[ E_{21} = \frac{28 \times 30}{60} = 14 \]
  • Untuk ASI Non Eksklusif dan Gizi Baik: \[ E_{22} = \frac{28 \times 30}{60} = 14 \]
  1. Chi-Square Statistic

Statistik Chi-Square dihitung menggunakan rumus: \[ \chi^2 = \sum \frac{(O - E)^2}{E} \] Di mana \(O\) adalah nilai observasi dan \(E\) adalah nilai yang diharapkan.

  • Untuk ASI Eksklusif dan Gizi Kurang: \[ \frac{(9 - 16)^2}{16} = \frac{49}{16} \approx 3.0625 \]
  • Untuk ASI Eksklusif dan Gizi Baik: \[ \frac{(23 - 16)^2}{16} = \frac{49}{16} \approx 3.0625 \]
  • Untuk ASI Non Eksklusif dan Gizi Kurang: \[ \frac{(21 - 14)^2}{14} = \frac{49}{14} \approx 3.5 \]
  • Untuk ASI Non Eksklusif dan Gizi Baik: \[ \frac{(7 - 14)^2}{14} = \frac{49}{14} \approx 3.5 \] Sehingga, statistik Chi-Square-nya adalah: \[ \chi^2 \approx 3.0625 + 3.0625 + 3.5 + 3.5 = 13.125 \]

Perhitungan R :

# Membuat tabel kontingensi
observed <- matrix(c(9, 23, 21, 7), nrow = 2, byrow = TRUE)
dimnames(observed) <- list("ASI" = c("Eksklusif", "Non Eksklusif"),
                           "Gizi" = c("Kurang", "Baik"))
cat("Observed Counts:\n")
## Observed Counts:
print(observed)
##                Gizi
## ASI             Kurang Baik
##   Eksklusif          9   23
##   Non Eksklusif     21    7
# Menghitung total baris, kolom, dan grand total
row_totals <- rowSums(observed)
col_totals <- colSums(observed)
grand_total <- sum(observed)

# Menghitung expected counts secara manual
expected <- outer(row_totals, col_totals, FUN = function(r, c) r * c / grand_total)
cat("\nExpected Counts:\n")
## 
## Expected Counts:
print(expected)
##               Kurang Baik
## Eksklusif         16   16
## Non Eksklusif     14   14
# Menghitung Chi-Square secara manual
chi_sq_manual <- sum((observed - expected)^2 / expected)
cat("\nManual Chi-square Statistic:", round(chi_sq_manual, 4), "\n")
## 
## Manual Chi-square Statistic: 13.125
# Menggunakan fungsi chisq.test untuk verifikasi
chi_test <- chisq.test(observed, correct = FALSE)  # non-corrected for continuity
cat("\nOutput dari chisq.test:\n")
## 
## Output dari chisq.test:
print(chi_test)
## 
##  Pearson's Chi-squared test
## 
## data:  observed
## X-squared = 13.125, df = 1, p-value = 0.0002914

Dari hasil perhitungan, didapatkan hasil nilai Chi - square yaitu 13,125 dan nilai p - value sebesar 0,0002. Jika p -value dibandingkan dengan taraf signifikansi 0,05 maka keputusan yang diambil adalah tolak H0 . Artinya, terdapat hubungan antara variabel “Pemberian ASI” dengan “Kondisi Gizi”.

6.2.4 Uji Exact Fisher

Uji Exact Fisher digunakan untuk menguji hubungan antara dua variabel kategorikal dalam tabel kontingensi kecil, dimana asumsi Chi - square tidak berlaku karena ukuran sampel yang kecil. Metode ini dapat menjadi allternatif uji Chi - Square terutama pada ukuran sampel yang kecil, Uji Exact Fisher tidak memerlukan asumsi normalitas serta memberikan hasil dengan keakuratan yang lebih tinggi dibandingkan menggunakan Chi - Square dengan jumlah frekuensi yang lebih kecil.

  1. Perhitungan Manual

Pada uji Fisher, probabilitas untuk mendapatkan tabel teramati dihitung dengan distribusi hipergeometrik. Untuk tabel 2x2, probabilitas tabel teramati adalah:

\[ P = \frac{\binom{32}{9} \, \binom{28}{21}}{\binom{60}{30}} \]

di mana: - \(\displaystyle \binom{32}{9}\) adalah jumlah cara memilih 9 kejadian (gizi kurang) dari 32 individu pada kelompok ASI Eksklusif. - \(\displaystyle \binom{28}{21}\) adalah jumlah cara memilih 21 kejadian (gizi kurang) dari 28 individu pada kelompok ASI Non Eksklusif. - \(\displaystyle \binom{60}{30}\) adalah total cara memilih 30 kasus gizi kurang dari 60 individu.

Perhitungan manual p-value Fisher secara lengkap memerlukan penjumlahan semua probabilitas tabel yang sama atau lebih ekstrem dari tabel yang teramati. Di bawah ini, sebagai ilustrasi, kita hitung probabilitas untuk tabel yang teramati :

# Menghitung probabilitas tabel teramati secara manual:
p_obs <- (choose(32, 9) * choose(28, 21)) / choose(60, 30)
cat("Probabilitas tabel teramati (P_obs):", p_obs, "\n")
## Probabilitas tabel teramati (P_obs): 0.0002808187
  1. Perhitungan dengan Package R:
# Membuat tabel kontingensi
table_data <- matrix(c(9, 23, 21, 7), nrow = 2, byrow = TRUE)
dimnames(table_data) <- list("ASI" = c("Eksklusif", "Non Eksklusif"),
                             "Gizi" = c("Kurang", "Baik"))

# Uji Fisher
fisher_result <- fisher.test(table_data)
print(fisher_result)
## 
##  Fisher's Exact Test for Count Data
## 
## data:  table_data
## p-value = 0.0006346
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
##  0.03487475 0.46857232
## sample estimates:
## odds ratio 
##  0.1357577

Didapatkan p - value sebesar 0,00063 < taraf signifikansi (0,05) yang artinya hasil ini sudah cukup untuk membuktikan adanya hubungan yang signifikan antara variabel “Pemberian ASI” dan “Kondisi Gizi”.

6.3 Analisis Residual dalam Tabel Kontingensi

Residual pada tabel kontingensi digunakan untuk mengidentifikasi sel-sel yang memiliki perbedaan signifikan antara frekuensi yang diobservasi dan frekuensi yang diharapkan; dengan kata lain, residual menunjukkan selisih antara nilai aktual dan nilai yang diperhitungkan berdasarkan model independensi. Jika nilai residual suatu sel jauh dari nol, entah itu positif ataupun negatif, ini menandakan bahwa frekuensi observasi sangat menyimpang dari ekspektasi, sehingga sel tersebut memberikan kontribusi yang besar terhadap asosiasi antar variabel. Sebaliknya, jika nilai residual mendekati nol, ini menunjukkan bahwa nilai observasi hampir sama dengan nilai yang diharapkan, sehingga sel tersebut tidak banyak berperan dalam hubungan antar variabel.

Contohnya:

  • Residual ≈ 0: Menunjukkan tidak adanya perbedaan signifikan antara nilai yang diobservasi dan yang diharapkan. Hal ini diartikan bahwa variabel baris dan kolom tidak menunjukkan keterkaitan yang kuat.

  • Residual positif besar: Frekuensi observasi jauh lebih tinggi daripada yang diharapkan, mencerminkan adanya hubungan positif yang signifikan.

  • Residual negatif besar: Frekuensi observasi jauh lebih rendah dari nilai ekspektasi, mengindikasikan hubungan negatif atau kurangnya asosiasi.

Jenis Residual

Dalam analisis tabel kontingensi terdapat dua jenis residual yang umum digunakan, yaitu:

  1. Pearson Residual

Pearson Residual mengukur deviasi antara nilai observasi (\(O_{ij}\)) dan nilai yang diharapkan (\(E_{ij}\)) relatif terhadap akar variansi yang diasumsikan sama dengan nilai yang diharapkan (tanpa penyesuaian tambahan). Rumusnya adalah:

\[ r_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}}} \]

Dimana: - \(O_{ij}\) = frekuensi yang diobservasi pada sel ke-\(ij\). - \(E_{ij}\) = frekuensi yang diharapkan pada sel ke-\(ij\), yang dihitung dengan:

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

Pearson Residual memberikan informasi dasar mengenai seberapa jauh penyimpangannya dari nilai yang diharapkan dalam satuan standar deviasi.

Perhitungan R:

# Data observasi
observed <- matrix(c(9, 23, 21, 7), nrow = 2, byrow = TRUE)
rownames(observed) <- c("Ekslusif", "Non Ekslusif")
colnames(observed) <- c("Gizi Kurang", "Gizi Baik")

# Hitung total baris, kolom, dan grand total
row_totals <- rowSums(observed)
col_totals <- colSums(observed)
grand_total <- sum(observed)

# Hitung expected counts
expected <- outer(row_totals, col_totals, FUN = function(r, c) r * c / grand_total)

# Hitung Pearson Residual
pearson_residual <- (observed - expected) / sqrt(expected)

cat("Observed Counts:\n")
## Observed Counts:
print(observed)
##              Gizi Kurang Gizi Baik
## Ekslusif               9        23
## Non Ekslusif          21         7
cat("\nExpected Counts:\n")
## 
## Expected Counts:
print(expected)
##              Gizi Kurang Gizi Baik
## Ekslusif              16        16
## Non Ekslusif          14        14
cat("\nPearson Residuals:\n")
## 
## Pearson Residuals:
print(pearson_residual)
##              Gizi Kurang Gizi Baik
## Ekslusif       -1.750000  1.750000
## Non Ekslusif    1.870829 -1.870829
  1. Standardized Residual

    Standardized Residual memberikan ukuran yang lebih terstandarisasi dengan memperhitungkan proporsi total baris dan kolom dalam tabel. Rumus yang sering digunakan adalah:

\[ r^*_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}\left(1 - \frac{\text{Total Baris}_i}{\text{Grand Total}}\right)\left(1 - \frac{\text{Total Kolom}_j}{\text{Grand Total}}\right)}} \]

Dimana: - \(\frac{\text{Total Baris}_i}{\text{Grand Total}}\) adalah proporsi total pada baris ke-\(i\). - \(\frac{\text{Total Kolom}_j}{\text{Grand Total}}\) adalah proporsi total pada kolom ke-\(j\).

Standardized Residual membantu mengidentifikasi sel-sel yang secara relatif memiliki kontribusi yang lebih besar terhadap ketidaksesuaian model independensi, sehingga memudahkan interpretasi apakah suatu sel sangat menyimpang atau tidak.

Perhitungan R:

# Membuat tabel kontingensi dengan data ASI
observed <- matrix(c(9, 23, 21, 7), nrow = 2, byrow = TRUE)
rownames(observed) <- c("ASI Eksklusif", "ASI Non Eksklusif")
colnames(observed) <- c("Gizi Kurang", "Gizi Baik")

# Hitung total baris, total kolom, dan grand total
row_totals <- rowSums(observed)
col_totals <- colSums(observed)
grand_total <- sum(observed)

# Hitung expected counts
expected <- outer(row_totals, col_totals, FUN = function(r, c) r * c / grand_total)

# Inisialisasi matriks untuk standardized residual
std_residual <- matrix(NA, nrow = nrow(observed), ncol = ncol(observed))

# Menghitung standardized residual untuk tiap sel
for (i in 1:nrow(observed)) {
  for (j in 1:ncol(observed)) {
    se_adjusted <- sqrt(expected[i, j] * (1 - (row_totals[i] / grand_total)) * (1 - (col_totals[j] / grand_total)))
    std_residual[i, j] <- (observed[i, j] - expected[i, j]) / se_adjusted
  }
}

# Memberikan nama pada baris dan kolom matriks residual
rownames(std_residual) <- rownames(observed)
colnames(std_residual) <- colnames(observed)

# Menampilkan hasil
cat("Tabel Observasi:\n")
## Tabel Observasi:
print(observed)
##                   Gizi Kurang Gizi Baik
## ASI Eksklusif               9        23
## ASI Non Eksklusif          21         7
cat("\nTabel Expected Counts:\n")
## 
## Tabel Expected Counts:
print(expected)
##                   Gizi Kurang Gizi Baik
## ASI Eksklusif              16        16
## ASI Non Eksklusif          14        14
cat("\nStandardized Residuals:\n")
## 
## Standardized Residuals:
print(round(std_residual, 2))
##                   Gizi Kurang Gizi Baik
## ASI Eksklusif           -3.62      3.62
## ASI Non Eksklusif        3.62     -3.62

7 Tabel Kontingensi 3 Arah

Dalam praktik penelitian, tabel kontingensi tiga arah tak sekadar menyajikan angka, melainkan mengisahkan interaksi dinamis antara variabel-variabel kategori. Misalnya, dalam kajian epidemiologi, penelitian sering mengeksplorasi hubungan antara kebiasaan merokok (variabel 𝑋) dan kejadian kanker paru-paru (variabel 𝑌) dengan mengendalikan peran usia (variabel 𝑍). Dengan demikian, peneliti dapat mengisolasi efek murni merokok terhadap kesehatan paru-paru—menyingkirkan kerancuan yang mungkin muncul akibat perbedaan usia. Analisis kondisional seperti ini membuka jendela pemahaman yang lebih jernih mengenai bagaimana faktor risiko benar-benar bekerja.

Sementara itu, dalam dunia penelitian sosial dan hukum, tabel kontingensi tiga arah juga digunakan untuk menyelidiki bagaimana variabel-variabel saling memengaruhi. Contohnya, dalam menilai dampak ras tersangka (variabel 𝑋) terhadap keputusan hukuman mati (variabel 𝑌), peneliti dapat mempertimbangkan peran ras korban (variabel 𝑍) sebagai faktor penentu. Pendekatan ini tidak hanya mengungkap pola secara agregat, tetapi juga memaparkan lapisan kompleks bias dan pertimbangan etis yang tersembunyi, sehingga memberikan wawasan mendalam yang tak akan terlihat bila hanya mengandalkan data marginal.

Lebih jauh lagi, inovasi dalam narasi analisis seperti ini juga mendorong aplikasinya ke bidang lain misalnya pemasaran. Di sini, hubungan antara preferensi konsumen (variabel 𝑋) dan keputusan pembelian (variabel 𝑌) seringkali dipengaruhi oleh demografi seperti usia atau pendapatan (variabel 𝑍). Dengan menggunakan tabel parsial—yang menyuguhkan analisis kondisional—peneliti dapat menggali insight strategis yang lebih presisi, mengantisipasi pergeseran tren serta kebutuhan pasar yang spesifik di setiap segmen.

Pendekatan inovatif ini, dengan menekankan peran variabel pengganggu, membantu menghindari distorsi interpretasi seperti yang disorot oleh Simpson’s Paradox. Sehingga, dibandingkan dengan melihat data secara agregat melalui tabel marginal, analisis menggunakan tabel parsial memberikan resolusi lebih tinggi atas kerumitan hubungan antar variabel.

# Membuat array 3 dimensi untuk tabel kontingensi 3 arah
data_array <- array(
  c(
    # X2 = sangat kurang
    3, 24, 7,      # Y = cukup, baik, sangat baik | X1 = sangat kurang
    2, 18, 16,     # Y = cukup, baik, sangat baik | X1 = kurang
    0, 11, 18,     # Y = cukup, baik, sangat baik | X1 = cukup
    0, 1, 0,       # Y = cukup, baik, sangat baik | X1 = baik
    # X2 = kurang
    3, 37, 25,
    2, 10, 9,
    0, 7, 7,
    0, 0, 0,       # X1 = baik tidak ada data untuk X2 = kurang
    # X2 = cukup
    0, 7, 7,
    0, 0, 0,
    0, 0, 0,
    0, 0, 0
  ),
  dim = c(3, 4, 3),
  dimnames = list(
    Y = c("cukup", "baik", "sangat baik"),
    X1 = c("sangat kurang", "kurang", "cukup", "baik"),
    X2 = c("sangat kurang", "kurang", "cukup")
  )
)

# Menampilkan tabel
data_array
## , , X2 = sangat kurang
## 
##              X1
## Y             sangat kurang kurang cukup baik
##   cukup                   3      2     0    0
##   baik                   24     18    11    1
##   sangat baik             7     16    18    0
## 
## , , X2 = kurang
## 
##              X1
## Y             sangat kurang kurang cukup baik
##   cukup                   3      2     0    0
##   baik                   37     10     7    0
##   sangat baik            25      9     7    0
## 
## , , X2 = cukup
## 
##              X1
## Y             sangat kurang kurang cukup baik
##   cukup                   0      0     0    0
##   baik                    7      0     0    0
##   sangat baik             7      0     0    0

7.1 Tabel Parsial dan Marginal

7.1.1 Tabel Parsial

Tabel parsial adalah tabel kontingensi dua arah (dalam kasus ini 3×3) yang memperlihatkan distribusi frekuensi antara variabel X dan Y untuk masing-masing kategori tertentu dari variabel Z. Dengan kata lain, setiap tabel parsial mengisolasi data untuk satu tingkat Z saja sehingga kita bisa melihat bagaimana hubungan X dan Y berlangsung ketika faktor Z dikendalikan.

Misalkan variabel Z memiliki dua kategori, katakanlah Z₁ dan Z₂. Maka, kita memiliki dua tabel parsial:

Tabel Parsial untuk Z = Z₁

Tabel Parsial untuk Z = Z₂

Dalam masing-masing tabel, misalnya nilai n₁₂¹ menunjukkan jumlah kasus di mana X = X₁ dan Y = Y₂ ketika Z bernilai Z₁.

Tabel parsial memungkinkan analisis hubungan bersyarat antara X dan Y. Dengan mengontrol Z (misalnya, usia, ras, atau karakteristik lain), peneliti dapat mengevaluasi apakah asosiasi antara X dan Y stabil atau berubah-ubah di setiap kategori Z. Cara ini membantu mengidentifikasi fenomena seperti Simpson’s Paradox, di mana pola asosiasi agregat bisa berbeda dengan pola pada tingkat subkelompok.

7.1.2 Tabel Marginal

Tabel marginal merupakan tabel 3×3 yang diperoleh dengan mengakumulasi (menjumlahkan) data dari semua kategori pada variabel Z. Dengan mengabaikan atau “menghapus” detail tentang Z, kita melihat gambaran keseluruhan (agregat) dari hubungan antara X dan Y.

Tabel marginal memberikan gambaran agregat tentang hubungan antara X dan Y tanpa memisahkan efek dari variabel Z. Meskipun ini berguna untuk memperoleh intuisi global, pengabaian Z seringkali menyembunyikan variasi atau bahkan mengaburkan hubungan sebenarnya antara X dan Y. Itulah sebabnya, dalam analisis yang lebih mendalam, terutama ketika potensi pengaruh faktor perancu ada, tabel parsial biasanya menjadi fokus utama.

Implementasi R:

# Membuat data frame kontingensi 3 arah
data_kontingensi <- expand.grid(
  Z = c("Perkotaan", "Pedesaan"),  # Perkotaan, Pedesaan
  X = c("rendah", "menengah", "tinggi"), #Tingkat pendidikan
  Y = c("tidak bekerja", "bekerja informal", "bekerja formal") #Pekerjaan
)

# Menambahkan frekuensi sesuai dengan tabel di atas
data_kontingensi$Freq <- c(
  # Z = Perkotaan 
  10, 25, 5,
  6, 15, 20,
  1, 4, 30,
  # Z = Pedesaaan
  20, 30, 2,
  12, 18, 8,
  4, 6, 10
)

# Menampilkan dataset
print(data_kontingensi)
##            Z        X                Y Freq
## 1  Perkotaan   rendah    tidak bekerja   10
## 2   Pedesaan   rendah    tidak bekerja   25
## 3  Perkotaan menengah    tidak bekerja    5
## 4   Pedesaan menengah    tidak bekerja    6
## 5  Perkotaan   tinggi    tidak bekerja   15
## 6   Pedesaan   tinggi    tidak bekerja   20
## 7  Perkotaan   rendah bekerja informal    1
## 8   Pedesaan   rendah bekerja informal    4
## 9  Perkotaan menengah bekerja informal   30
## 10  Pedesaan menengah bekerja informal   20
## 11 Perkotaan   tinggi bekerja informal   30
## 12  Pedesaan   tinggi bekerja informal    2
## 13 Perkotaan   rendah   bekerja formal   12
## 14  Pedesaan   rendah   bekerja formal   18
## 15 Perkotaan menengah   bekerja formal    8
## 16  Pedesaan menengah   bekerja formal    4
## 17 Perkotaan   tinggi   bekerja formal    6
## 18  Pedesaan   tinggi   bekerja formal   10

7.2 Distribusi Peluang

7.2.1 Peluang Bersama

Untuk mencari nilai peluang bersama digunakan rumus :

\[ P(X = x, Y = y) = \frac{\sum_{z} n_{xyz}}{N} \]

Dengan :

  • \(n_{xyz}\): Frekuensi observasi untuk kombinasi \(X = x\), \(Y = y\), dan \(Z = z\).
  • \(N\): Total frekuensi semua observasi pada tabel kontingensi.

Perhitungan R:

# Membuat data frame kontingensi 3 arah
data_kontingensi <- expand.grid(
  Z = c("Perkotaan", "Pedesaan"),
  X = c("rendah", "menengah", "tinggi"),
  Y = c("tidak bekerja", "bekerja informal", "bekerja formal")
)

# Menambahkan frekuensi sesuai dengan tabel
data_kontingensi$Freq <- c(
  # Untuk Z = Perkotaan 
  10, 25, 5,     # X = rendah, Y = "tidak bekerja", "bekerja informal", "bekerja formal"
   6, 15, 20,    # X = menengah, ..
   1, 4, 30,     # X = tinggi, ..
  # Untuk Z = Pedesaan
  20, 30, 2,    # X = rendah, ..
  12, 18, 8,    # X = menengah, ..
   4,  6, 10    # X = tinggi, ..
)

# Menampilkan dataset
print(data_kontingensi)
##            Z        X                Y Freq
## 1  Perkotaan   rendah    tidak bekerja   10
## 2   Pedesaan   rendah    tidak bekerja   25
## 3  Perkotaan menengah    tidak bekerja    5
## 4   Pedesaan menengah    tidak bekerja    6
## 5  Perkotaan   tinggi    tidak bekerja   15
## 6   Pedesaan   tinggi    tidak bekerja   20
## 7  Perkotaan   rendah bekerja informal    1
## 8   Pedesaan   rendah bekerja informal    4
## 9  Perkotaan menengah bekerja informal   30
## 10  Pedesaan menengah bekerja informal   20
## 11 Perkotaan   tinggi bekerja informal   30
## 12  Pedesaan   tinggi bekerja informal    2
## 13 Perkotaan   rendah   bekerja formal   12
## 14  Pedesaan   rendah   bekerja formal   18
## 15 Perkotaan menengah   bekerja formal    8
## 16  Pedesaan menengah   bekerja formal    4
## 17 Perkotaan   tinggi   bekerja formal    6
## 18  Pedesaan   tinggi   bekerja formal   10
# Total frekuensi, yang akan digunakan sebagai penyebut peluang
total_N <- sum(data_kontingensi$Freq)
cat("Total frekuensi (N):", total_N, "\n")
## Total frekuensi (N): 226
# Mendapatkan frekuensi untuk sel X = "rendah", Y = "tidak bekerja", Z = "Perkotaan"
n_xyz <- data_kontingensi$Freq[data_kontingensi$X == "rendah" & 
                                data_kontingensi$Y == "tidak bekerja" & 
                                data_kontingensi$Z == "Perkotaan"]

# Peluang gabungan
p_xyz_manual <- n_xyz / total_N
cat("P(X = rendah, Y = tidak bekerja, Z = Perkotaan) =", p_xyz_manual, "\n")
## P(X = rendah, Y = tidak bekerja, Z = Perkotaan) = 0.04424779

7.2.2 Peluang Marginal

Peluang marginal dihitung dengan rumus :

\[ P(X = x, Y = y) = \frac{\sum_{z} n_{xyz}}{N} \]

Perhitungan R:

# Menghitung total frekuensi untuk X = "rendah" dan Y = "tidak bekerja" dari semua kategori Z
n_xy <- sum(data_kontingensi$Freq[data_kontingensi$X == "rendah" & 
                                   data_kontingensi$Y == "tidak bekerja"])
p_xy_manual <- n_xy / total_N
cat("P(X = rendah, Y = tidak bekerja) =", p_xy_manual, "\n")
## P(X = rendah, Y = tidak bekerja) = 0.1548673

7.2.3 Peluang Bersyarat

Peluang bersyarat dapat dicari apabila peluang suatu kejadian mensyaratkan suatu kejadian terjadi lebih dulu. Dihitung dengan rumus :

\[ P(Y = y \mid X = x, Z = z) = \frac{n_{xyz}}{\sum_{y} n_{xyz}} \]

Perhitungan R:

# Numerator: frekuensi untuk X = "rendah", Y = "tidak bekerja", Z = "Perkotaan"
n_numerator <- data_kontingensi$Freq[data_kontingensi$X == "rendah" &
                                      data_kontingensi$Y == "tidak bekerja" &
                                      data_kontingensi$Z == "Perkotaan"]

# Denominator: total frekuensi untuk X = "rendah" dan Z = "Perkotaan"
n_denominator <- sum(data_kontingensi$Freq[data_kontingensi$X == "rendah" &
                                            data_kontingensi$Z == "Perkotaan"])

p_cond_manual <- n_numerator / n_denominator
cat("P(Y = tidak bekerja | X = rendah, Z = Perkotaan) =", p_cond_manual, "\n")
## P(Y = tidak bekerja | X = rendah, Z = Perkotaan) = 0.4347826
# Menggunakan tabel parsial untuk Z = "Perkotaan"
partial_perkotaan <- xtabs(Freq ~ X + Y, data = subset(data_kontingensi, Z == "Perkotaan"))
p_cond_R <- partial_perkotaan["rendah", "tidak bekerja"] / sum(partial_perkotaan["rendah",])
cat("P(Y = tidak bekerja | X = rendah, Z = Perkotaan) =", p_cond_R, "\n")
## P(Y = tidak bekerja | X = rendah, Z = Perkotaan) = 0.4347826

7.3 Ukuran Asosiasi

  • Risk Difference (RD): Selisih risiko antara kelompok eksposur dan non-eksposur, yaitu

\[ RD = \text{Risk}_{\text{eksposur}} - \text{Risk}_{\text{non-eksposur}} \]

  • Risk Ratio (RR): Perbandingan risiko antara kelompok eksposur dan non-eksposur, yaitu

\[ RR = \frac{\text{Risk}_{\text{eksposur}}}{\text{Risk}_{\text{non-eksposur}}} \]

  • Odds Ratio (OR): Perbandingan odds (peluang) kejadian antara kelompok eksposur dan non-eksposur, yaitu

\[ OR = \frac{\text{Odds}_{\text{eksposur}}}{\text{Odds}_{\text{non-eksposur}}} \]

Karena data awal merupakan tabel kontingensi tiga arah, akan dilakukan dikotomisasi data agar analisis 2×2 dapat dilakukan dengan:

  • Variabel Z: Berfokus pada wilayah Perkotaan.

  • Variabel X: Membandingkan kelompok “rendah” (eksposur) dengan “tinggi” (non-eksposur).

  • Variabel Y: Pendefinian Event jika “tidak bekerja” dan Non-Event jika “bekerja informal” atau “bekerja formal”.

Risk Difference (RD)

Risk pada suatu kelompok adalah:

\[ \text{Risk} = \frac{\text{Jumlah Event}}{\text{Total Observasi pada kelompok tersebut}} \]

Sehingga,

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

di mana: - \(a\) = jumlah event pada kelompok eksposur (rendah)
- \(b\) = jumlah non-event pada kelompok eksposur
- \(c\) = jumlah event pada kelompok non-eksposur (tinggi)
- \(d\) = jumlah non-event pada kelompok non-eksposur

Risk Ratio (RR)

Menunjukkan seberapa besar risiko kejadian “tidak bekerja” pada kelompok “rendah” dibandingkan dengan “tinggi”

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

Odds Ratio (OR)

Mengukur perbandingan peluang terjadinya “tidak bekerja” antara kedua kelompok. Odds dihitung sebagai:

\[ \text{Odds} = \frac{\text{Jumlah Event}}{\text{Jumlah Non-Event}} \]

Sehingga,

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

Perhitungan R:

# Membuat data frame kontingensi 3 arah
data_kontingensi <- expand.grid(
  Z = c("Perkotaan", "Pedesaan"),
  X = c("rendah", "menengah", "tinggi"),
  Y = c("tidak bekerja", "bekerja informal", "bekerja formal")
)

# Menambahkan frekuensi sesuai dengan tabel yang diberikan
data_kontingensi$Freq <- c(
  # Data untuk Z = Perkotaan 
  10, 25, 5,    # X = rendah, masing-masing Y
   6, 15, 20,   # X = menengah, masing-masing Y
   1, 4, 30,    # X = tinggi, masing-masing Y
  # Data untuk Z = Pedesaan
  20, 30, 2,    # X = rendah
  12, 18, 8,    # X = menengah
   4,  6, 10    # X = tinggi
)

# Tampilkan data
data_kontingensi
# Subset data untuk wilayah Perkotaan dan X "rendah" atau "tinggi"
data_sub <- subset(data_kontingensi, Z == "Perkotaan" & X %in% c("rendah", "tinggi"))

# Recode variabel Y menjadi biner: "Event" untuk "tidak bekerja", "Non-Event" untuk sisanya
data_sub$Y_binary <- ifelse(data_sub$Y == "tidak bekerja", "Event", "Non-Event")

# Agregasi data untuk membuat tabel 2×2
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.3.3
## 
## 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
library(tidyr)
## Warning: package 'tidyr' was built under R version 4.3.3
table_2x2 <- data_sub %>%
  group_by(X, Y_binary) %>%
  summarise(Freq = sum(Freq)) %>%
  pivot_wider(names_from = Y_binary, values_from = Freq)
## `summarise()` has grouped output by 'X'. You can override using the `.groups`
## argument.
table_2x2
# Konversi tabel ke data frame agar mudah mengakses nilai
table_df <- as.data.frame(table_2x2)

# Pisahkan data untuk kelompok eksposur (rendah) dan non-eksposur (tinggi)
exposed <- table_df[table_df$X == "rendah", ]
non_exposed <- table_df[table_df$X == "tinggi", ]

# Misalkan:
# a = Event pada kelompok "rendah"
# b = Non-Event pada kelompok "rendah"
a <- exposed$Event
b <- exposed$`Non-Event`
n_exposed <- a + b

# c = Event pada kelompok "tinggi"
# d = Non-Event pada kelompok "tinggi"
c <- non_exposed$Event
d <- non_exposed$`Non-Event`
n_non_exposed <- c + d

# Hitung Risiko untuk masing-masing kelompok
risk_exposed <- a / n_exposed
risk_non_exposed <- c / n_non_exposed

cat("Risk pada kelompok 'rendah':", risk_exposed, "\n")
## Risk pada kelompok 'rendah': 0.4347826
cat("Risk pada kelompok 'tinggi':", risk_non_exposed, "\n")
## Risk pada kelompok 'tinggi': 0.2941176
# Risk Difference (RD)
RD <- risk_exposed - risk_non_exposed

# Risk Ratio (RR)
RR <- risk_exposed / risk_non_exposed

# Odds untuk masing-masing kelompok
odds_exposed <- a / b
odds_non_exposed <- c / d

# Odds Ratio (OR)
OR <- odds_exposed / odds_non_exposed

cat("Risk Difference (RD):", RD, "\n")
## Risk Difference (RD): 0.140665
cat("Risk Ratio (RR):", RR, "\n")
## Risk Ratio (RR): 1.478261
cat("Odds Ratio (OR):", OR, "\n")
## Odds Ratio (OR): 1.846154

7.3.1 Tabel Kontingensi Parsial

Tabel kontingensi parsial digunakan untuk melihat hubungan antara dua variabel (misalnya, \(X\) dan \(Y\)) pada tiap-tiap level variabel ketiga (\(Z\)).
Artinya, untuk setiap nilai \(z\) pada variabel \(Z\), dapat dibuat matriks frekuensi \(X \times Y\) yang merepresentasikan distribusi data kondisi \(Z=z\).

Misalkan \(n_{xyz}\) adalah frekuensi observasi pada kombinasi \(X=x\), \(Y=y\), dan \(Z=z\). Maka tabel parsial untuk level \(Z=z\) adalah:

\[ T^{(z)}_{xy} = n_{xy}(z) = n_{xyz} \quad \text{dengan total}\quad N(z)=\sum_{x}\sum_{y} n_{xyz}. \]

Dapat dihitung distribusi proporsi relatif dengan:

\[ P(X=x, Y=y \mid Z=z) = \frac{n_{xyz}}{N(z)}. \]

Dengan menggunakan data kontingensi 3 arah dengan tiga variabel: - Z: “Perkotaan” dan “Pedesaan” - X: “rendah”, “menengah”, “tinggi” - Y: “tidak bekerja”, “bekerja informal”, “bekerja formal”

Implementasi R:

# Membuat data frame kontingensi 3 arah
data_kontingensi <- expand.grid(
  Z = c("Perkotaan", "Pedesaan"),
  X = c("rendah", "menengah", "tinggi"),
  Y = c("tidak bekerja", "bekerja informal", "bekerja formal")
)

# Menambahkan frekuensi berdasarkan data
data_kontingensi$Freq <- c(
  # Data untuk Z = Perkotaan 
  10, 25, 5,    # X = rendah dengan masing-masing kategori Y
   6, 15, 20,   # X = menengah
   1, 4, 30,    # X = tinggi
  # Data untuk Z = Pedesaan
  20, 30, 2,    # X = rendah 
  12, 18, 8,    # X = menengah
   4,  6, 10    # X = tinggi
)

# Tampilkan data
data_kontingensi
# Subset data untuk Z = "Perkotaan"
data_perkotaan <- subset(data_kontingensi, Z == "Perkotaan")

# Membuat tabel kontingensi parsial antara variabel X dan Y untuk wilayah Perkotaan
tabel_parsial <- xtabs(Freq ~ X + Y, data = data_perkotaan)
# Perbaikan: pastikan nama variabel konsisten, gunakan data_perkotaan
tabel_parsial <- xtabs(Freq ~ X + Y, data = data_perkotaan)

# Tampilkan tabel parsial
tabel_parsial
##           Y
## X          tidak bekerja bekerja informal bekerja formal
##   rendah              10                1             12
##   menengah             5               30              8
##   tinggi              15               30              6
# Menghitung total frekuensi pada masing-masing tabel parsial (untuk Z = Perkotaan)
total_perkotaan <- sum(tabel_parsial)
# Menghitung proporsi relatif
proporsi_parsial <- tabel_parsial / total_perkotaan

proporsi_parsial
##           Y
## X          tidak bekerja bekerja informal bekerja formal
##   rendah     0.085470085      0.008547009    0.102564103
##   menengah   0.042735043      0.256410256    0.068376068
##   tinggi     0.128205128      0.256410256    0.051282051

7.4 Conditional Independence

Dalam analisis data kategorikal yang melibatkan tiga variabel – misalnya, \(X\), \(Y\), dan \(Z\)conditional independence atau independensi bersyarat antara \(X\) dan \(Y\) diberikan \(Z\) menyatakan bahwa informasi mengenai \(X\) tidak memberikan tambahan informasi mengenai \(Y\) apabila kita sudah mengetahui \(Z\). Secara matematis, \(X\) dan \(Y\) dikatakan independen bersyarat pada \(Z\) jika, untuk setiap nilai \(x\), \(y\), dan \(z\),

\[ P(X=x, Y=y \mid Z=z) = P(X=x \mid Z=z) \times P(Y=y \mid Z=z). \]

Jika rumus di atas terpenuhi pada setiap level \(z\), maka tidak ada asosiasi langsung antara \(X\) dan \(Y\) ketika efek \(Z\) telah dikontrol.

Studi Kasus: Tabel Kontingensi 3 Arah

Misalkan terdapati data dengan tiga variabel: - \(Z\): Lokasi (contoh: “Perkotaan” dan “Pedesaan”) - \(X\): Kategori, misalnya tingkat pendidikan (“rendah”, “menengah”, “tinggi”) - \(Y\): Status pekerjaan (“tidak bekerja”, “bekerja informal”, “bekerja formal”)

Untuk lebih jelas, akan dilakukan perhitungan di R:

# Subset data untuk Z = "Perkotaan"
data_urban <- subset(data_kontingensi, Z == "Perkotaan")

# Buat tabel kontingensi parsial antara X dan Y untuk Z = "Perkotaan"
table_urban <- xtabs(Freq ~ X + Y, data = data_urban)
table_urban
##           Y
## X          tidak bekerja bekerja informal bekerja formal
##   rendah              10                1             12
##   menengah             5               30              8
##   tinggi              15               30              6
# Hitung joint probability P(X, Y | Z = Perkotaan)
prop_joint <- prop.table(table_urban)
cat("Joint Probability (P(X, Y | Z = Perkotaan)):\n")
## Joint Probability (P(X, Y | Z = Perkotaan)):
print(prop_joint)
##           Y
## X          tidak bekerja bekerja informal bekerja formal
##   rendah     0.085470085      0.008547009    0.102564103
##   menengah   0.042735043      0.256410256    0.068376068
##   tinggi     0.128205128      0.256410256    0.051282051
# Hitung probabilitas marginal untuk X dan Y masing-masing, diberikan Z = Perkotaan
row_margin <- margin.table(table_urban, margin = 1)  # Total untuk setiap X
col_margin <- margin.table(table_urban, margin = 2)  # Total untuk setiap Y
total_urban <- sum(table_urban)

prop_X <- row_margin / total_urban
prop_Y <- col_margin / total_urban

# Jika X dan Y independen bersyarat, maka:
# P(X=x, Y=y | Z=Perkotaan) = P(X=x | Z=Perkotaan) * P(Y=y | Z=Perkotaan)
prod_matrix <- outer(prop_X, prop_Y, FUN = "*")
cat("\nProduk Probabilitas Marginal (P(X|Z) * P(Y|Z)):\n")
## 
## Produk Probabilitas Marginal (P(X|Z) * P(Y|Z)):
print(prod_matrix)
##           Y
## X          tidak bekerja bekerja informal bekerja formal
##   rendah      0.05040544        0.1024911     0.04368471
##   menengah    0.09423625        0.1916137     0.08167142
##   tinggi      0.11176857        0.2272628     0.09686610

Pengujian independensi juga dapat dilakukan dengan statistik Uji Chi - Square

# Menambahkan frekuensi berdasarkan data contoh
data_kontingensi$Freq <- c(
  # Untuk Z = Perkotaan
  10, 25, 5,    # X = rendah dengan masing-masing Y
   6, 15, 20,   # X = menengah
   1,  4, 30,   # X = tinggi
  # Untuk Z = Pedesaan
  20, 30, 2,    # X = rendah
  12, 18, 8,    # X = menengah
   4,  6, 10    # X = tinggi
)

# Tampilkan data
data_kontingensi
# Membuat tabel marginal dengan mengakumulasi semua kategori Z
table_marginal <- xtabs(Freq ~ X + Y, data = data_kontingensi)
print(table_marginal)
##           Y
## X          tidak bekerja bekerja informal bekerja formal
##   rendah              35                5             30
##   menengah            11               50             12
##   tinggi              35               32             16
# Melakukan uji chi-square pada tabel marginal
chi_test_marginal <- chisq.test(table_marginal)
print(chi_test_marginal)
## 
##  Pearson's Chi-squared test
## 
## data:  table_marginal
## X-squared = 60.373, df = 4, p-value = 2.421e-12

Hasil uji chi-square pada setiap tabel parsial menguji apakah terdapat asosiasi antara

𝑋 dan 𝑌 dalam konteks spesifik 𝑍. Dengan membandingkan hasil untuk masing-masing level

𝑍(misalnya, Perkotaan vs. Pedesaan), kita dapat mengevaluasi apakah hubungan antara 𝑋 dan 𝑌 konsisten di seluruh tingkat 𝑍 atau ada perbedaan signifikan. Setelah dilakukan pengujian menggunakan statistik Uji Chi Square maka didapatkan hasil p - value yang signifikan, artinya “Status pekerjaan” dan “Tingkat pendidikan” masih berhubungan satu sama lain walaupun variabel Z yaitu lokasi tempat tinggal dikendalikan.

7.6 Inferensi Tabel Kontingensi 3 Arah

Tabel kontingensi tiga arah digunakan untuk menganalisis hubungan antara dua variabel kategorik dengan mempertimbangkan variabel kontrol. Sebagai contoh, pengkajian hubungan antara kebiasaan merokok (\(X\)) dan kanker paru-paru (\(Y\)) dengan mengendalikan pengaruh usia sebagai variabel \(Z\).

Data disusun sedemikian rupa sehingga untuk setiap kategori usia (\(Z\)) terdapat sebuah tabel parsial (2×2) yang menunjukkan distribusi frekuensi antara \(X\) dan \(Y\). Tabel marginal yang mengabaikan \(Z\) juga dapat digunakan, namun informasi dari variabel kontrol dapat menjadi kunci dalam menginterpretasikan asosiasi antara \(X\) dan \(Y\). Untuk setiap tabel parsial (2×2), odds ratio (OR) dihitung dengan rumus:

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

di mana: - \(a\) = jumlah kasus dimana \(X=1\) dan \(Y=1\)
- \(b\) = jumlah kasus dimana \(X=1\) dan \(Y=0\)
- \(c\) = jumlah kasus dimana \(X=0\) dan \(Y=1\)
- \(d\) = jumlah kasus dimana \(X=0\) dan \(Y=0\)

Jika odds ratio parsial, yakni OR yang dihitung untuk tiap level \(Z\), relatif konstan di seluruh strata, kita dapat menghitung odds ratio bersama menggunakan estimasi Mantel-Haenszel dengan rumus:

\[ OR_{MH} = \frac{\sum_{z} \frac{a_z \, d_z}{N_z}}{\sum_{z} \frac{b_z \, c_z}{N_z}} \]

di mana: - \(a_z\), \(b_z\), \(c_z\), \(d_z\) adalah frekuensi pada masing-masing sel tabel untuk strata \(z\) (misalnya, kelompok usia tertentu), - \(N_z\) adalah total observasi pada strata \(z\).

Rumus tersebut memberikan ukuran asosiasi yang disesuaikan dengan variabel kontrol usia, sehingga memberikan gambaran yang lebih valid tentang hubungan antara \(X\) dan \(Y\) setelah mengendalikan perbedaan distribusi usia.

Misalkan terdapat data dengan tiga variabel: - \(X\): Kebiasaan merokok (1: Merokok, 0: Tidak Merokok) - \(Y\): Kanker paru-paru (1: Ada, 0: Tidak Ada) - \(Z\): Kelompok usia (misalnya, “Muda”, “Paruh Baya”, “Tua”)

Untuk setiap level \(Z\) kita dapat membuat tabel 2×2 yang berisikan jumlah frekuensi sesuai kombinasi \(X\) dan \(Y\). Perhitungan odds ratio pada setiap level \(Z\) dilakukan menggunakan rumus di atas. Jika nilai OR pada setiap strata cukup homogen, kita dapat menghitung \(OR_{MH}\) sebagai ukuran effect bersama.

Kesimpulan

Pendekatan ini memungkinkan peneliti:

  • Mengontrol variabel usia (\(Z\)): sehingga asosiasi \(X\) (merokok) dan \(Y\) (kanker paru-paru) dievaluasi dalam konteks yang lebih homogen.

  • Menghitung odds ratio parsial untuk masing-masing tingkat usia dan kemudian menggabungkannya dengan estimasi Mantel-Haenszel untuk mendapatkan \(OR_{MH}\).

  • Menafsirkan asosiasi: Jika OR parsial relatif konstan di semua strata, \(OR_{MH}\) adalah estimasi yang baik dari ukuran efek bersama \(X\) terhadap \(Y\) setelah mengendalikan usia.

7.6.1 Independensi Bersyarat dalam Tabel Kontingensi Tiga Arah

Independensi Bersyarat: Dua variabel \(X\) dan \(Y\) dikatakan independen bersyarat terhadap variabel \(Z\) jika, di setiap strata \(z\), rasio odds (odds ratio) antara \(X\) dan \(Y\) adalah 1. Secara matematis, kondisi ini dapat dituliskan sebagai:

\[ OR(X, Y \mid Z=z) = 1 \quad \text{untuk setiap } z. \]

Artinya, bila kita mengendalikan \(Z\), tidak terdapat asosiasi antara \(X\) dan \(Y\) dalam setiap kelompok \(z\).

\[ OR_z = \frac{a_z \times d_z}{b_z \times c_z}. \]

Dalam konteks independensi bersyarat, jika \(X\) dan \(Y\) benar-benar tidak berhubungan setelah mengontrol \(Z\), maka diharapkan:

\[ OR_z = 1 \quad \text{(untuk setiap } z\text{)}. \]

Selain itu, uji Cochran-Mantel-Haenszel (CMH) dapat digunakan untuk mengestimasi odds ratio gabungan di antara berbagai strata \(Z\) dengan rumus:

\[ OR_{MH} = \frac{\sum_{z} \frac{a_z d_z}{N_z}}{\sum_{z} \frac{b_z c_z}{N_z}}, \]

di mana \(N_z = a_z + b_z + c_z + d_z\) merupakan total pengamatan pada strata \(z\).

7.6.2 Pengujian Statistik untuk Independensi Bersyarat

Metode Cochran-Mantel-Haenszel (CMH) adalah teknik statistik yang digunakan untuk menguji asosiasi antara dua variabel biner yang diukur pada sejumlah strata (tingkat) yang ditentukan oleh variabel ketiga. Dengan kata lain, CMH menggabungkan informasi dari beberapa tabel kontingensi 2×2 (parsial) yang diperoleh dari masing-masing kelompok (misalnya, berdasarkan lokasi, kategori usia, atau faktor lain) untuk mendapatkan estimasi odds ratio (OR) yang disesuaikan dengan variabel perancu.

Dalam konteks data pekerjaan, misalkan ingin dilakukan pengkajian hubungan antara tingkat pendidikan (\(X\)) dan status pekerjaan (\(Y\)) dengan mempertimbangkan lokasi (\(Z\)). Di sini, kita akan membandingkan dua kelompok dalam \(X\) – misalnya “rendah” dan “tinggi” – serta mengklasifikasikan status pekerjaan \(Y\) menjadi dua kategori:
- Tidak Bekerja (jika nilai \(Y\) adalah “tidak bekerja”), dan
- Bekerja (menggabungkan “bekerja informal” dan “bekerja formal”).

CMH memungkinkan kita menguji hipotesis bahwa hubungan antara \(X\) dan \(Y\) tidak ada setelah dikontrol oleh \(Z\).

Hipotesis dan Statistik Uji

Hipotesis Nol (H₀):
Odds ratio (OR) antara \(X\) dan \(Y\) adalah sama dengan 1 pada setiap strata \(Z\); yaitu, tidak ada asosiasi antara \(X\) dan \(Y\) setelah mengendalikan \(Z\).

\[ H_0: OR(X, Y\mid Z=z) = 1 \quad untuk\ setiap\ z. \]

Hipotesis Alternatif (H₁):
Terdapat asosiasi antara \(X\) dan \(Y\) dalam setidaknya salah satu strata \(Z\) sehingga OR tidak sama dengan 1.

Statistik Uji:
Metode CMH menghitung sebuah statistik uji yang mengikuti distribusi chi-square dengan 1 derajat kebebasan (df). Nilai p yang signifikan (misalnya, p < 0.05) menolak hipotesis nol, yang berarti terdapat asosiasi yang konsisten antara \(X\) dan \(Y\) setelah mengendalikan \(Z\).

Selain itu, odds ratio gabungan dapat dihitung menggunakan rumus Mantel-Haenszel:

\[ OR_{MH} = \frac{\sum_{z} \dfrac{a_z \, d_z}{N_z}}{\sum_{z} \dfrac{b_z \, c_z}{N_z}}, \]

di mana:
- \(a_z\) = jumlah individu dengan \(X=1\) dan \(Y=1\) pada strata \(z\),
- \(b_z\) = jumlah individu dengan \(X=1\) dan \(Y=0\) pada strata \(z\),
- \(c_z\) = jumlah individu dengan \(X=0\) dan \(Y=1\) pada strata \(z\),
- \(d_z\) = jumlah individu dengan \(X=0\) dan \(Y=0\) pada strata \(z\),
- \(N_z = a_z+b_z+c_z+d_z\) adalah total observasi pada strata \(z\).

Perhitungan R:

# Data frame kontingensi 3 arah
data_kontingensi <- expand.grid(
  Z = c("Perkotaan", "Pedesaan"),
  X = c("rendah", "menengah", "tinggi"),
  Y = c("tidak bekerja", "bekerja informal", "bekerja formal")
)

# Menambahkan frekuensi (data contoh)
data_kontingensi$Freq <- c(
  # Data untuk Z = Perkotaan 
  10, 25, 5,    # X = rendah
   6, 15, 20,   # X = menengah
   1,  4, 30,   # X = tinggi
  # Data untuk Z = Pedesaan
  20, 30, 2,    # X = rendah
  12, 18, 8,    # X = menengah
   4,  6, 10    # X = tinggi
)

# Kita hanya gunakan dua kategori pendidikan untuk analisis 2x2: "rendah" dan "tinggi"
data_pekerjaan <- subset(data_kontingensi, X %in% c("rendah", "tinggi"))

# Mendichotomisasi status pekerjaan: 
# "tidak bekerja" tetap, "bekerja informal" & "bekerja formal" digabungkan menjadi "Bekerja"
data_pekerjaan$Status <- ifelse(data_pekerjaan$Y == "tidak bekerja", 
                                "Tidak Bekerja", "Bekerja")

# Tampilan data yang sudah didichotomisasi
data_pekerjaan
library(dplyr)
library(tidyr)

# Agregasi data agar memperoleh tabel 2x2 untuk setiap level Z
table_pekerjaan <- data_pekerjaan %>%
  group_by(Z, X, Status) %>%
  summarise(Freq = sum(Freq)) %>%
  pivot_wider(names_from = Status, values_from = Freq)
## `summarise()` has grouped output by 'Z', 'X'. You can override using the
## `.groups` argument.
print(table_pekerjaan)
## # A tibble: 4 × 4
## # Groups:   Z, X [4]
##   Z         X      Bekerja `Tidak Bekerja`
##   <fct>     <fct>    <dbl>           <dbl>
## 1 Perkotaan rendah      13              10
## 2 Perkotaan tinggi      36              15
## 3 Pedesaan  rendah      22              25
## 4 Pedesaan  tinggi      12              20
# Membuat tabel kontingensi 2x2 untuk setiap strata menggunakan xtabs
tabel_array <- xtabs(Freq ~ X + Status + Z, data = data_pekerjaan)
print(tabel_array)
## , , Z = Perkotaan
## 
##           Status
## X          Bekerja Tidak Bekerja
##   rendah        13            10
##   menengah       0             0
##   tinggi        36            15
## 
## , , Z = Pedesaan
## 
##           Status
## X          Bekerja Tidak Bekerja
##   rendah        22            25
##   menengah       0             0
##   tinggi        12            20

Kesimpulan Dari hasil uji Cochran-Mantel-Haenszel, kita dapat menilai apakah terdapat hubungan signifikan antara tingkat pendidikan dan status pekerjaan setelah mengontrol faktor Lokasi tempat tinggal. Jika nilai p lebih kecil dari 0.05, maka hipotesis nol ditolak, yang berarti ada hubungan signifikan antara status pekerjaan dan tingkat pendidikan setelah mempertimbangkan faktor lokasi

8 Generalized Linear Model (GLM)

Generalized Linear Model (GLM) merupakan perluasan dari regresi linear konvensional yang memungkinkan pemodelan variabel dependen yang tidak memenuhi asumsi normalitas. GLM mengintegrasikan tiga komponen utama:

1. Komponen Acak (Random Component): Mengasumsikan bahwa variabel respons \(Y\) mengikuti distribusi dari kelas exponential family (misalnya, binomial, Poisson, Gamma, dll).

2. Komponen Sistematis (Systematic Component): Membangun prediktor linier, yaitu kombinasi linear dari variabel bebas, yang dituliskan sebagai

\[ \eta = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \cdots + \beta_p X_p. \]

3. Fungsi Link (Link Function): Memetakan nilai harapan \(\mu = E(Y)\) ke prediktor linier \(\eta\) dengan fungsi link \(g(\cdot)\) sehingga

\[ g(\mu) = \eta. \]

8.1 Exponential Family

Distribusi dalam keluarga eksponensial memiliki bentuk umum:

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

di mana: - \(\theta\) adalah parameter natural (natural parameter), - \(\phi\) adalah parameter dispersi, - \(a(\phi)\) adalah fungsi skala, - \(b(\theta)\) berhubungan dengan fungsi cumulant yang menentukan nilai harapan dan variansi, - \(c(y,\phi)\) adalah fungsi normalisasi yang memastikan densitas terintegrasi menjadi 1.

Rumus inilah yang mendasari asumsi distribusi pada GLM. Misalnya, untuk distribusi binomial, fungsi-fungsi tersebut memiliki bentuk khusus sehingga menghasilkan link logit; sedangkan untuk distribusi Poisson, link log sering digunakan.

8.2 Regresi Logistik

Regresi logistik adalah metode statistik yang digunakan untuk memodelkan hubungan antara satu atau beberapa variabel independen dengan variabel dependen biner (dua kategori). Model ini menghubungkan nilai probabilitas kejadian sebuah peristiwa dengan prediktor melalui fungsi logit, yang merupakan transformasi dari fungsi sigmoid.

Rumus dasar model regresi logistik untuk satu prediktor adalah:

\[ \text{logit}(p) = \log\left(\frac{p}{1-p}\right) = \beta_0 + \beta_1 X, \]

atau secara eksplisit untuk probabilitas:

\[ p(X) = \frac{1}{1 + \exp\left[-(\beta_0 + \beta_1 X)\right]}. \]

Asumsi Regresi Logistik

Beberapa asumsi utama dalam regresi logistik meliputi:

  • Outcome Biner: Variabel dependen harus memiliki dua kategori (misalnya, 1 dan 0).
  • Observasi Independen: Data yang digunakan diharapkan tidak memiliki ketergantungan antar observasi.
  • Hubungan Linear pada Skala Logit: Logit dari probabilitas diharapkan memiliki hubungan linear dengan prediktor.
  • Tidak Ada Multikolinearitas: Variabel-variabel independen tidak seharusnya sangat berkorelasi satu sama lain.
  • Ukuran Sampel yang Memadai: Untuk menghasilkan estimasi parameter yang stabil, harus ada jumlah observasi yang cukup dalam setiap kategori outcome.

Estimasi Parameter

Parameter \(\beta_0\) dan \(\beta_1\) diestimasi menggunakan metode maksimum likelihood. Proses ini mencari nilai parameter yang memaksimalkan fungsi likelihood, yaitu fungsi yang mengukur seberapa “mungkin” data yang diamati diberikan parameter yang diusulkan.

Simulasi :

Untuk memahami bagaimana fungsi sigmoid bekerja, kita dapat mensimulasikan grafik fungsi sigmoid. Simulasi ini memperlihatkan bagaimana probabilitas berubah seiring dengan peningkatan nilai \(X\).

# Mendefinisikan fungsi sigmoid
sigmoid <- function(x, beta0 = 0, beta1 = 1) {
  1 / (1 + exp(-(beta0 + beta1 * x)))
}

# Membuat vektor nilai x
x_vals <- seq(-10, 10, length.out = 200)

# Menghitung probabilitas menggunakan berbagai kombinasi parameter
p_default <- sigmoid(x_vals)
p_shifted <- sigmoid(x_vals, beta0 = -2, beta1 = 1.5)

# Plot fungsi sigmoid
plot(x_vals, p_default, type = "l", lwd = 2, col = "blue",
     ylim = c(0, 1), xlab = "X", ylab = "Probabilitas p(X)",
     main = "Simulasi Fungsi Sigmoid")
lines(x_vals, p_shifted, lwd = 2, col = "red")
legend("topleft", legend = c("Default: β₀=0, β₁=1", "Shifted: β₀=-2, β₁=1.5"),
       col = c("blue", "red"), lwd = 2)

Sebagai contoh, kita akan menggunakan data mtcars untuk memodelkan probabilitas sebuah mobil memiliki transmisi manual (variabel am) berdasarkan variabel mpg dan wt. Pada contoh ini, variabel am (0 = otomatis, 1 = manual) merupakan outcome biner.

Perhitungan R

# Menggunakan dataset mtcars
data(mtcars)

# Ubah variabel 'am' menjadi faktor
mtcars$am <- factor(mtcars$am, labels = c("Otomatis", "Manual"))

# Mendefinisikan model regresi logistik
model_logistic <- glm(am ~ mpg + wt, data = mtcars, family = binomial)

# Menampilkan ringkasan model
summary(model_logistic)
## 
## Call:
## glm(formula = am ~ mpg + wt, family = binomial, data = mtcars)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)  
## (Intercept)  25.8866    12.1935   2.123   0.0338 *
## mpg          -0.3242     0.2395  -1.354   0.1759  
## wt           -6.4162     2.5466  -2.519   0.0118 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 43.230  on 31  degrees of freedom
## Residual deviance: 17.184  on 29  degrees of freedom
## AIC: 23.184
## 
## Number of Fisher Scoring iterations: 7

Ringkasan model menunjukkan estimasi parameter 𝛽 beserta nilai p yang dapat digunakan untuk menguji signifikansi efek prediktor terhadap probabilitas outcome.

8.3 Regresi Poisson

Regresi Poisson adalah jenis Generalized Linear Model (GLM) yang digunakan untuk memodelkan variabel dependen berupa count data (data hitungan), seperti jumlah kejadian suatu peristiwa dalam periode waktu atau ruang tertentu. Variabel response \(Y\) diasumsikan mengikuti distribusi Poisson, yang memiliki properti bahwa rata-rata dan variansnya sama.

Model Poisson menggunakan fungsi link logaritma, sehingga hubungan antara variabel dependen dan prediktor dapat dituliskan sebagai:

\[ \log(\mu) = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \cdots + \beta_p X_p, \]

atau secara implisit menyatakan bahwa:

\[ \mu = \exp(\beta_0 + \beta_1 X_1 + \cdots + \beta_p X_p), \]

di mana \(\mu = E(Y)\) merupakan rata-rata kejadian.

Asumsi Regresi Poisson

Beberapa asumsi penting pada regresi Poisson adalah:

  • Data Count: Variabel dependen berupa hitungan (non-negatif) seperti 0, 1, 2, …
  • Distribusi Poisson: \(Y\) diasumsikan mengikuti distribusi Poisson, sehingga \(Var(Y) = E(Y) = \mu\).
  • Fungsi Link Logaritma: Fungsi link yang digunakan adalah log, sehingga hubungan antara prediktor dan nilai harapan bersifat multiplikatif.
  • Independensi Observasi: Observasi harus independen satu sama lain.
  • Tidak Ada Overdispersi: Varians seharusnya mendekati mean. Jika varians jauh lebih besar (overdispersi), model alternatif (seperti model negative binomial) mungkin lebih tepat.

Estimasi Parameter

Parameter \(\beta_0, \beta_1, \ldots, \beta_p\) diestimasi menggunakan metode maksimum likelihood. Pendekatan ini mencari nilai parameter yang memaksimalkan fungsi likelihood, mehngukur seberapa “mungkin” data yang diamati diberikan model yang diusulkan.

Simulasi Fungsi Poisson

Untuk memahami bagaimana perubahan nilai prediktor \(X\) mempengaruhi rata-rata kejadian, kita dapat mensimulasikan fungsi Poisson. Di bawah ini, kita mensimulasikan data count \(Y\) dengan asumsi model:

\[ \mu = \exp(\beta_0 + \beta_1 X), \]

misalnya dengan \(\beta_0 = 1\) dan \(\beta_1 = 0.5\).

Simulasi R:

set.seed(123)

# Simulasikan data
n <- 200
x <- runif(n, 0, 10)  # nilai X dari distribusi uniform
beta0 <- 1
beta1 <- 0.5
# Menghitung lambda untuk setiap nilai X
lambda_sim <- exp(beta0 + beta1 * x)
# Mensimulasikan data count Y dari distribusi Poisson
y <- rpois(n, lambda_sim)

# Gabungkan data ke dalam data frame
sim_data <- data.frame(x = x, y = y)

# Plot data simulasi
plot(sim_data$x, sim_data$y, xlab = "X", ylab = "Count Y",
     main = "Data Simulasi Count", col = "blue", pch = 16)

# Fit model regresi Poisson
model_pois <- glm(y ~ x, data = sim_data, family = poisson)
summary(model_pois)
## 
## Call:
## glm(formula = y ~ x, family = poisson, data = sim_data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) 0.947438   0.035665   26.57   <2e-16 ***
## x           0.507015   0.004329  117.13   <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: 20785.11  on 199  degrees of freedom
## Residual deviance:   208.11  on 198  degrees of freedom
## AIC: 1276.1
## 
## Number of Fisher Scoring iterations: 4

Pada ringkasan model, kita dapat melihat estimasi parameter untuk intercept ( 𝛽0) dan slope (𝛽1), nilai p, dan ukuran goodness-of-fit. Parameter-parameter tersebut diinterpretasikan dalam skala log, sehingga kenaikan satu unit di 𝑋 berhubungan dengan perubahan exp⁡(𝛽1)pada rata-rata count 𝑌 .

9 Inferensi GLM

9.1 Ekspektasi dan Varians dalam GLM

Distribusi dari keluarga eksponensial dapat dituliskan secara umum sebagai

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

di mana: - \(\theta\) adalah parameter natural, - \(\phi\) adalah parameter dispersi, - \(a(\phi)\) adalah fungsi skala, - \(b(\theta)\) adalah fungsi cumulant, - \(c(y,\phi)\) adalah fungsi normalisasi.

Dari bentuk tersebut, kita dapat menurunkan:

  • Ekspektasi:

\[ E(Y) = b'(\theta) \]

  • Varians:

\[ Var(Y) = b''(\theta) \, a(\phi) \]

Sebagai contoh, untuk regresi logistik (keluarga binomial) maupun regresi Poisson, fungsi \(a(\phi)\) dan \(b(\theta)\) memiliki bentuk khusus sehingga menghasilkan fungsi link dan bentuk varians yang sesuai.

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 tidak sama dengan konstan

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

Contoh:

• Poisson: 𝑉(𝜇) = 𝜇

• Binomial: 𝑉(𝜇) = 𝜇(1−𝜇)

Contoh Regresi Poisson

# Simulasi data
set.seed(123)
x <- rnorm(100)
mu <- exp(0.5 + 0.8 * x)
y <- rpois(100, mu)
model <- glm(y ~ x, family = poisson)
summary(model)
## 
## Call:
## glm(formula = y ~ x, family = poisson)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.44950    0.08872   5.066 4.05e-07 ***
## x            0.86000    0.07463  11.523  < 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: 245.05  on 99  degrees of freedom
## Residual deviance: 106.78  on 98  degrees of freedom
## AIC: 325.76
## 
## Number of Fisher Scoring iterations: 5
  • Estimasi Maksimum Likelihood (MLE)
    Prinsip utama: memaksimalkan fungsi likelihood atau log-likelihood.
    Tahapan:
    Susun fungsi likelihood \(L(\theta | x)\) berdasarkan distribusi probabilitas data sampel, di mana \(\theta\) adalah parameter yang akan diestimasi dan \(x\) adalah data sampel.

  • Ambil logaritma natural dari fungsi likelihood untuk mendapatkan log-likelihood, \(\ell(\theta) = \ln L(\theta | x)\), untuk mempermudah perhitungan.

  • Hitung turunan pertama log-likelihood terhadap parameter \(\theta\), yakni \(\frac{d\ell(\theta)}{d\theta}\).

  • Setel turunan pertama sama dengan nol, \(\frac{d\ell(\theta)}{d\theta} = 0\), untuk menemukan nilai \(\theta\) yang memaksimalkan fungsi.

  • Gunakan turunan kedua, \(\frac{d^2\ell(\theta)}{d\theta^2}\), untuk memastikan bahwa solusi adalah maksimum (harus negatif di titik solusi).

  • Substitusikan solusi \(\theta\) ke dalam persamaan untuk mendapatkan estimasi maksimum likelihood.


Misalkan \(x_1, x_2, \dots, x_n\) adalah sampel acak dari distribusi dengan fungsi densitas probabilitas \(f(x | \theta)\). Maka:

Fungsi likelihood:

\[ L(\theta | x) = \prod_{i=1}^n f(x_i | \theta) \]

Fungsi log-likelihood:

\[ \ell(\theta) = \ln L(\theta | x) = \sum_{i=1}^n \ln f(x_i | \theta) \]

Estimasi MLE diperoleh dengan memecahkan:

\[ \frac{d\ell(\theta)}{d\theta} = 0 \]

Karena model GLM tidak memiliki bentuk eksplisit, metode numerik digunakan untuk penyelesaian.

Metode Optimisasi Newton-Raphson
Prinsip : Memanfaatkan vektor skor (gradien) & Menggunakan matriks Hessian untuk optimisasi.

Tahapan :

  • Pilih fungsi tujuan \(f(\theta)\) yang akan dioptimalkan (biasanya log-likelihood untuk estimasi).

  • Tentukan turunan pertama (gradien) \(\nabla f(\theta) = \frac{df(\theta)}{d\theta}\) dan turunan kedua (Hessian) \(\frac{d^2f(\theta)}{d\theta^2}\).

  • Tentukan nilai awal parameter \(\theta_0\).

  • Gunakan rumus iterasi Newton-Raphson: \[ \theta_{n+1} = \theta_n - \left( \frac{d^2f(\theta_n)}{d\theta^2} \right)^{-1} \frac{df(\theta_n)}{d\theta} \]

  • Ulangi langkah 4 hingga konvergensi, yaitu ketika \(|\theta_{n+1} - \theta_n| < \epsilon\) (dengan \(\epsilon\) kecil) atau gradien mendekati nol.

  • Pastikan solusi adalah maksimum dengan memeriksa bahwa matriks Hessian negatif definit (untuk maksimasi) di titik solusi.

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.

Kesimpulan:

• Ekspektasi digunakan untuk mengetahui ketakbiasan estimasi

• Varians digunakan untuk mengukur presisi dan menyusun uji statistik

• Distribusi asimptotik dari𝛽 sangat bergantung pada kedua konsep ini

• Dalam GLM, varians sangat tergantung pada bentuk distribusi eksponensial dari data

9.2 Penaksiran Parameter dalam GLM

Metode Penaksiran Parameter: Maximum Likelihood Estimation (MLE)

Dalam GLM, parameter \(\beta\) diestimasi dengan memaksimalkan fungsi log-likelihood

\[ \ell(\beta) = \sum_{i=1}^{n} \log f_Y(y_i; \theta_i, \phi), \]

dengan hubungan \(\theta_i\) yang berkaitan dengan \(X_i\) melalui fungsi link. Proses optimasi umumnya dilakukan dengan algoritma iteratif (misalnya, iteratively reweighted least squares, IRLS).

Metode Penaksiran Parameter : Newton-Raphson

Metode Newton-Raphson dimanfaatkan untuk menentukan nilai parameter \(\beta\) yang memaksimalkan fungsi log-likelihood dalam model regresi logistik.

Untuk menjelaskan langkah-langkah estimasi parameter menggunakan metode Newton-Raphson dalam konteks memaksimalkan fungsi log-likelihood (misalnya, pada regresi Poisson atau regresi logistik) dan menyusunnya dalam LaTeX, berikut adalah penjelasan dan kode LaTeX yang sesuai.

Langkah-langkah Metode Newton-Raphson untuk Estimasi Parameter:

1. Inisialisasi: Tentukan nilai awal untuk vektor parameter, misalnya \(\beta^{(0)}\).

2. Hitung Vektor Skor dan Matriks Hessian:

  • Vektor skor \(U(\beta^{(t)})\) adalah gradien (turunan pertama) dari fungsi log-likelihood terhadap \(\beta\).

  • Matriks Hessian \(H(\beta^{(t)})\) adalah matriks turunan kedua dari fungsi log-likelihood terhadap \(\beta\).

3. Perbarui Parameter: Gunakan rumus iterasi Newton-Raphson: \[ \beta^{(t+1)} = \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}). \] 4. Cek Konvergensi: Ulangi langkah 2 dan 3 hingga perubahan pada \(\beta\) cukup kecil (misalnya, \(|\beta^{(t+1)} - \beta^{(t)}| < \epsilon\)) atau log-likelihood konvergen.

5. Hasil Akhir: Nilai \(\beta^{(t+1)}\) yang konvergen adalah estimasi MLE untuk parameter \(\beta\).

set.seed(1)
x <- rnorm(100)
beta_true <- c(-1, 2)
X <- cbind(1, x)
eta <- X %*% beta_true
p <- 1 / (1 + exp(-eta))
y <- rbinom(100, 1, p)

9.3 Diagnostik Model GLM

Diagnostik Model GLM

Setelah model diestimasikan, diagnostik model memegang peranan evaluasi kelayakan model. Beberapa metode diagnostik meliputi:

  • Residual Analysis:
    Menghitung residual seperti deviance residual, Pearson residual, dan standardized residual untuk mengidentifikasi outlier atau ketidaksesuaian model.
  • Goodness-of-Fit:
    Menggunakan deviance atau Pearson chi-square statistic untuk menilai seberapa baik model cocok dengan data.
    Contoh:

\[ \text{Deviance} = 2 \sum_{i=1}^{n} \left[\ell(y_i; y_i) - \ell(y_i; \hat{\mu}_i)\right]. \]

  • Pengaruh Pengamatan:
    Identifikasi leverage, Cook’s distance, dan pengaruh individual observasi terhadap parameter model.

Ekspektasi dan Varians Estimator

Dengan menggunakan MLE, kita memperoleh estimator parameter \(\hat{\beta}\) yang, dengan asumsi model telah memenuhi syarat, memiliki properti:

  • Ekspektasi Estimator:
    Secara asimptotik, estimator MLE adalah tidak bias, yaitu:

\[ E(\hat{\beta}) \approx \beta. \]

  • Varians Estimator:
    Varians dari estimator MLE mendekati invers dari matriks Fisher:

\[ Var(\hat{\beta}) \approx I(\hat{\beta})^{-1}, \]

di mana \(I(\hat{\beta})\) merupakan matriks informasi Fisher, yaitu matriks kedua derivatif log-likelihood.

9.4 Metode Estimasi dan Inferensi Regresi Logistik

Regresi Logistik digunakan untuk outcome biner. Model logistik dinyatakan sebagai:

\[ \log\left(\frac{p}{1-p}\right) = \beta_0 + \beta_1 X_1 + \cdots + \beta_p X_p, \]

atau secara eksplisit untuk probabilitas

\[ p(X) = \frac{1}{1 + \exp\left[-(\beta_0 + \beta_1 X_1 + \cdots + \beta_p X_p)\right]}. \]

Estimasi Parameter:
Parameter diestimasi menggunakan MLE. Inferensi dilakukan menggunakan uji statistik (misalnya, uji Wald, uji likelihood ratio) dan interval kepercayaan diperoleh dari matriks informasi Fisher.

9.5 Metode Estimasi dan Inferensi Regresi Poisson

Regresi Poisson digunakan untuk memodelkan data hitungan (count data) di mana variabel dependen \(Y\) diasumsikan mengikuti distribusi Poisson. Model Poisson mengasumsikan bahwa rata-rata dan varians dari \(Y\) sama dengan \(\lambda\), yaitu:

\[ \log(\lambda_i) = \beta_0 + \beta_1 X_{i}, \]

atau secara eksplisit:

\[ \lambda_i = \exp\left(\beta_0 + \beta_1 X_i\right). \]

Dengan fungsi distribusi poisson :

\[P(X = k) = \frac{\lambda^k e^{-\lambda}}{k!}\]

Parameter \(\beta_0\) dan \(\beta_1\) diestimasi menggunakan Maximum Likelihood Estimation (MLE) melalui proses optimasi yang memaksimalkan log-likelihood fungsi dari model Poisson.

Simulasi Data

set.seed(123)
n <- 100
x <- rnorm(n)
X <- cbind(1, x) # Tambah intercept
beta_true <- c(0.5, 0.8)
eta <- X %*% beta_true
lambda <- exp(eta)
y <- rpois(n, lambda)
data<-data.frame(x=x, y=y)

Metode Newton Raphson

# Data simulasi
set.seed(1)
x <- rnorm(100)
beta_true <- c(-1, 2)
X <- cbind(1, x)
eta <- X %*% beta_true
p <- 1 / (1 + exp(-eta))
y <- rbinom(100, 1, p)

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- 7
beta # hasil estimasi
##         [,1]
##   -1.6325651
## x  0.9552515

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)  -1.6326     0.2517  -6.486 8.80e-11 ***
## x             0.9553     0.2076   4.602 4.18e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 72.613  on 99  degrees of freedom
## Residual deviance: 50.901  on 98  degrees of freedom
## AIC: 116.9
## 
## Number of Fisher Scoring iterations: 5
  • IRLS menyediakan pendekatan iteratif untuk menghitung estimasi MLE pada regresi Poisson.
  • Hasil perhitungan manual IRLS sangat mirip dengan keluaran fungsi glm() di R.
  • Metode ini membantu memahami secara mendalam mekanisme yang mendasari fungsi glm().

Pengujian hipotesis Uji Wald

Untuk menguji H0: \(\beta_j\)= 0

Simulasi Uji Wald

set.seed(123)
n <- 100
x <- rnorm(n)
log_odds <--0.5 + 1.2 * x
p <- 1 / (1 + exp(-log_odds))
y <- rbinom(n, 1, p)
data <- data.frame(x, y)
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.3097     0.2296  -1.349    0.177    
## x             1.2663     0.3080   4.111 3.94e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 137.99  on 99  degrees of freedom
## Residual deviance: 114.76  on 98  degrees of freedom
## AIC: 118.76
## 
## Number of Fisher Scoring iterations: 4

Langkah 1: Ambil nilai koefisien dan 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 
## 4.110965

Langkah 3 : Hitung Statistik Uji Wald

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

Langkah 4 : Hitung p-value

p_value <- 1- pchisq(Wald_stat, df = 1)
p_value
##            x 
## 3.940095e-05

Interpretasi • Jika p-value < 0.05, maka koefisien signifikan → variabel prediktor berpengaruh. • Jika p-value > 0.05, maka tidak ada cukup bukti untuk menolak H0. 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 individualdalam model regresi logistik.

Uji Likelihood Ratio (Chi-Square)

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

Evaluasi Model (AIC & BIC)

AIC

AIC(model)
## [1] 118.7598

BIC

BIC(model)
## [1] 123.9701
  • Estimasi parameter dalam regresi Poisson dihitung dengan metode Maksimum Likelihood Estimation (MLE).

  • Pengujian hipotesis dilakukan menggunakan uji Wald dan uji Rasio Likelihood.

  • Kriteria AIC dan BIC digunakan untuk mengevaluasi dan memilih model yang paling optimal.

Estimasi Parameter dan Log-Likelihood

Fungsi likelihood untuk \(n\) observasi bernilai

\[ L(\beta) = \prod_{i=1}^{n} \frac{e^{-\lambda_i} \lambda_i^{y_i}}{y_i!}, \quad \text{dengan } \lambda_i = \exp\left(\beta_0 + \beta_1 X_i\right). \]

Log-likelihood yang harus dimaksimalkan adalah

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

Proses optimasi (biasanya menggunakan algoritma iteratively reweighted least squares, IRLS) menghasilkan estimator \(\hat{\beta}_0\) dan \(\hat{\beta}_1\) beserta standar error-nya.

Uji Wald

Uji Wald digunakan untuk menguji signifikansi masing-masing parameter. Untuk tiap parameter \(\beta_j\), statistik uji Wald didefinisikan sebagai

\[ z_j = \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)}, \]

dan nilainya dibandingkan dengan distribusi normal standar. Nilai p yang dihasilkan memberikan informasi apakah variabel \(X\) secara signifikan mempengaruhi response \(Y\) pada level yang ditetapkan.

Uji Likelihood Ratio (LRT)

Uji Likelihood Ratio membandingkan model penuh (dengan semua prediktor) dengan model yang dibatasi (misalnya model null yang hanya memiliki intercept). Statistika uji

\[ D = -2\left[\ell(\text{model reduce}) - \ell(\text{model penuh})\right], \]

mengikuti distribusi chi-square dengan derajat kebebasan yang sama dengan selisih jumlah parameter antara kedua model. Nilai p yang kecil (misalnya, \(p < 0.05\)) menolak hipotesis nol bahwa model yang lebih sederhana sudah memadai.

Simulasi R :

set.seed(123)
n <- 200
x <- runif(n, 0, 10)  # Variabel prediktor
beta0 <- 1           # Intercept
beta1 <- 0.5         # Koefisien prediktor

lambda <- exp(beta0 + beta1 * x)  # Rata-rata (λ)
y <- rpois(n, lambda)             # Data hitungan

# Gabungkan ke dalam data frame
sim_data <- data.frame(x = x, y = y)

# Plot data simulasi
plot(sim_data$x, sim_data$y, 
     xlab = "X", ylab = "Count Y",
     main = "Data Simulasi Count", col = "blue", pch = 16)

10 Regresi Logistik dengan Prediktor Nominal, Ordinal, dan Rasio

Bab ini akan mengkaji bagaimana regresi logistik digunakan sebagai pemodelan variabel dependen biner yang dipengaruhi oleh tiga tipe prediktor:

  • Nominal: Kategori tanpa urutan, seperti Jenis Kelamin (Laki-laki atau Perempuan).
  • Ordinal: Kategori berurutan tanpa jarak pasti, contohnya Tingkat Kepuasan (Rendah, Sedang, Tinggi).
  • Rasio: Variabel kontinu dengan nol mutlak dan selisih bermakna, misalnya Usia atau Pendapatan.

10.1 Simulasi Data

Dataset simulasi berisi 200 observasi dengan variabel:

  • lokasi_hunian: Nominal, dengan kategori urban, suburban, rural (rural sebagai referensi).
  • sikap (sikap kepedulian terhadap lingkungan): Ordinal, dengan tingkatan tidak_peduli, cukup_peduli, sangat_peduli.
  • pendapatan_tahunan: Rasio, dihasilkan dari distribusi normal (mean = 150 juta rupiah, sd = 40 juta rupiah).
  • keputusan_membeli: Biner, dihasilkan berdasarkan probabilitas dari kombinasi prediktor.
# Mengatur seed untuk reproduktifitas
set.seed(456)
# Membuat dataset simulasi
n <- 200
data <- tibble(
lokasi_hunian = sample(c("urban", "suburban", "rural"), n, replace = TRUE, prob = c(0.4, 0.35, 0.25)),
sikap = sample(c("tidak_peduli", "cukup_peduli", "sangat_peduli"), n, replace = TRUE, prob = c(0.3, 0.5, 0.2)),
pendapatan_tahunan = pmax(rnorm(n, mean = 150, sd = 40), 0),
keputusan_membeli = rbinom(n, 1, 0.3) # Probabilitas awal
)
# Menyesuaikan keputusan_membeli berdasarkan prediktor
data$sikap <- factor(data$sikap, levels = c("tidak_peduli", "cukup_peduli", "sangat_peduli"))
logit <- -5 +
  0.5 * (data$lokasi_hunian == "suburban") +
  1.0 * (data$lokasi_hunian == "urban") +
  0.8 * (as.numeric(data$sikap) - 1) +
  0.02 * data$pendapatan_tahunan
prob <- pmax(pmin(1 / (1 + exp(-logit)), 1), 0)
data$keputusan_membeli <- rbinom(n, 1, prob)
head(data)

10.2 Ekplorasi Data

data %>%
  group_by(keputusan_membeli) %>%
  summarise(Jumlah = n(),
            rata_rata_pendapatan = mean(pendapatan_tahunan))

Terlihat karakteristik dari kelompok “Membeli kendaraan listrik” dan “Tidak membeli kendaraan listrik”.

10.3 Perlakuan Variabel Ordinal

10.3.1 Treat sebagai Nominal

data_nominal <- data %>%
mutate(
sikap = factor(sikap, levels = c("tidak_peduli", "cukup_peduli", "sangat_peduli"))
)
model_nominal <- glm(keputusan_membeli ~ lokasi_hunian + sikap + pendapatan_tahunan, data = data_nominal, family = binomial)
summary(model_nominal)
## 
## Call:
## glm(formula = keputusan_membeli ~ lokasi_hunian + sikap + pendapatan_tahunan, 
##     family = binomial, data = data_nominal)
## 
## Coefficients:
##                        Estimate Std. Error z value Pr(>|z|)    
## (Intercept)           -4.253239   0.883423  -4.814 1.48e-06 ***
## lokasi_huniansuburban  0.296501   0.419528   0.707  0.47972    
## lokasi_hunianurban     1.028451   0.414230   2.483  0.01304 *  
## sikapcukup_peduli      0.828329   0.426311   1.943  0.05201 .  
## sikapsangat_peduli     1.919444   0.492258   3.899 9.65e-05 ***
## pendapatan_tahunan     0.014552   0.004509   3.227  0.00125 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 260.19  on 199  degrees of freedom
## Residual deviance: 224.97  on 194  degrees of freedom
## AIC: 236.97
## 
## Number of Fisher Scoring iterations: 4

Baseline yang Dipakai

Dalam model regresi logistik, baseline adalah kategori referensi untuk variabel kategorikal, yaitu kategori yang tidak muncul dalam output koefisien.

  • Untuk variabel lokasi_hunian, baseline adalah “rural” karena tidak ada koefisien untuk “rural” dalam output.

  • Untuk variabel sikap, baseline adalah “tidak_peduli” karena tidak ada koefisien untuk “tidak_peduli”.

Interpretasi Setiap Koefisien

Setiap koefisien dalam model regresi logistik merepresentasikan perubahan dalam log odds untuk membeli (yaitu, log dari rasio peluang membeli vs. tidak membeli) ketika variabel prediktor berubah, dengan variabel lain tetap.

(Intercept): -4.253239

  • Log Odds: Log odds untuk membeli bagi individu pada baseline (lokasi hunian “rural” dan sikap “tidak_peduli”) adalah -4.253239.
  • Odds Ratio: exp(-4.253239) ≈ 0.0142, yang berarti peluang membeli sangat rendah (1.42%) untuk baseline.
  • P-value: 1.48e-06 (sangat signifikan), menunjukkan bahwa intercept ini berbeda secara signifikan dari nol.

lokasi_huniansuburban: 0.296501

  • Log Odds: Dibandingkan dengan baseline “rural”, log odds untuk membeli bagi individu di lokasi “suburban” meningkat sebesar 0.296501.
  • Odds Ratio: exp(0.296501) ≈ 1.345, artinya peluang membeli 1.345 kali lebih tinggi di suburban dibandingkan rural.
  • P-value: 0.47972 (tidak signifikan pada α=0.05), menunjukkan bahwa perbedaan ini tidak cukup kuat untuk dianggap signifikan.

lokasi_hunianurban: 1.028451

  • Log Odds: Dibandingkan dengan baseline “rural”, log odds untuk membeli bagi individu di lokasi “urban” meningkat sebesar 1.028451.
  • Odds Ratio: exp(1.028451) ≈ 2.797, artinya peluang membeli hampir 2.8 kali lebih tinggi di urban dibandingkan rural.
  • P-value: 0.01304 (signifikan pada α=0.05), menunjukkan bahwa lokasi urban secara signifikan meningkatkan peluang membeli.

sikapcukup_peduli: 0.828329

  • Log Odds: Dibandingkan dengan baseline “tidak_peduli”, log odds untuk membeli bagi individu dengan sikap “cukup_peduli” meningkat sebesar 0.828329.
  • Odds Ratio: exp(0.828329) ≈ 2.289, artinya peluang membeli 2.289 kali lebih tinggi untuk sikap “cukup_peduli”.
  • P-value: 0.05201 (marginal signifikan pada α=0.05), menunjukkan bahwa pengaruhnya hampir signifikan, tetapi tidak sepenuhnya meyakinkan.

sikapsangat_peduli: 1.919444

  • Log Odds: Dibandingkan dengan baseline “tidak_peduli”, log odds untuk membeli bagi individu dengan sikap “sangat_peduli” meningkat sebesar 1.919444.
  • Odds Ratio: exp(1.919444) ≈ 6.817, artinya peluang membeli hampir 7 kali lebih tinggi untuk sikap “sangat_peduli”.
  • P-value: 9.65e-05 (sangat signifikan), menunjukkan bahwa sikap “sangat_peduli” secara kuat meningkatkan peluang membeli.

pendapatan_tahunan: 0.014552

  • Log Odds: Untuk setiap peningkatan satu unit (juta rupiah) dalam pendapatan tahunan, log odds untuk membeli meningkat sebesar 0.014552.
  • Odds Ratio: exp(0.014552) ≈ 1.015, artinya setiap peningkatan satu unit pendapatan tahunan (juta rupiah) meningkatkan peluang membeli sebesar 1.5%.
  • P-value: 0.00125 (sangat signifikan), menunjukkan bahwa pendapatan tahunan adalah prediktor yang kuat.

Signifikansi Model

Berdasarkan p-value dari masing-masing koefisien, kita dapat menentukan variabel mana yang signifikan dalam meningkatkan peluang membeli:

Variabel yang signifikan :

lokasi_hunianurban (p = 0.01304)

sikapsangat_peduli (p = 9.65e-05)

pendapatan_tahunan (p = 0.00125)

Variabel yang marginal signifikan (p ≈ 0.05):

sikapcukup_peduli (p = 0.05201), hampir signifikan.

Variabel yang tidak signifikan: - lokasi_huniansuburban (p = 0.47972)

Jadi, lokasi hunian urban, sikap sangat peduli, dan pendapatan tahunan secara signifikan meningkatkan peluang membeli, sedangkan sikap cukup peduli memiliki pengaruh yang hampir signifikan, dan lokasi suburban tidak berpengaruh signifikan.

Interpretasi Goodness-of-Fit

Goodness-of-fit menunjukkan seberapa baik model menjelaskan variasi dalam data:

Null Deviance: 260.19 (deviance dari model tanpa prediktor, hanya intercept).

Residual Deviance: 224.97 (deviance dari model dengan prediktor).

Penurunan Deviance: 260.19 - 224.97 = 35.22, yang menunjukkan bahwa model dengan prediktor lebih baik dalam menjelaskan data dibandingkan model null (tanpa prediktor).

AIC (Akaike Information Criterion): 236.97, digunakan untuk membandingkan model. Nilai AIC yang lebih rendah menunjukkan model yang lebih baik (akan dibandingkan dengan model selanjutnya).

Secara keseluruhan, model ini lebih baik daripada model null, tetapi residual deviance yang masih cukup besar (224.97) menunjukkan bahwa ada variasi dalam data yang belum dijelaskan oleh model ini. Ini berarti model mungkin dapat ditingkatkan dengan menambahkan prediktor lain atau mempertimbangkan interaksi antar variabel.

Kesimpulan Praktis

Berdasarkan hasil model regresi logistik:

Faktor yang Meningkatkan Peluang Membeli: -

  • Lokasi Hunian Urban: Individu di area urban memiliki peluang membeli yang signifikan lebih tinggi dibandingkan dengan area rural.

  • Sikap Sangat Peduli: Individu yang sangat peduli memiliki peluang membeli yang jauh lebih tinggi dibandingkan dengan yang tidak peduli.

  • Pendapatan Tahunan: Semakin tinggi pendapatan tahunan, semakin besar peluang untuk membeli.

  • Faktor Marginal: Sikap “cukup_peduli” hampir signifikan, menunjukkan bahwa individu dengan sikap ini mungkin juga lebih cenderung membeli, tetapi bukti tidak sekuat untuk sikap “sangat_peduli”.

Faktor Tidak Signifikan: Lokasi hunian suburban tidak menunjukkan perbedaan signifikan dalam peluang membeli dibandingkan dengan rural.

Secara keseluruhan, model ini memberikan wawasan yang berguna, tetapi ada ruang untuk perbaikan dengan mengeksplorasi variabel tambahan atau interaksi antar variabel untuk menjelaskan lebih banyak variasi dalam keputusan membeli.

10.3.2 Treat sebagai Numerik

L

data_numeric <- data %>%
mutate(
sikap_num = case_when(sikap == "tidak_peduli"  ~ 1,
sikap == "cukup_peduli"  ~ 2,
sikap == "sangat_peduli" ~ 3)
)
model_num <- glm(keputusan_membeli ~ lokasi_hunian + sikap_num + pendapatan_tahunan, data = data_numeric, family = binomial)
summary(model_num)
## 
## Call:
## glm(formula = keputusan_membeli ~ lokasi_hunian + sikap_num + 
##     pendapatan_tahunan, family = binomial, data = data_numeric)
## 
## Coefficients:
##                        Estimate Std. Error z value Pr(>|z|)    
## (Intercept)           -5.308233   0.985274  -5.388 7.14e-08 ***
## lokasi_huniansuburban  0.292033   0.418456   0.698  0.48525    
## lokasi_hunianurban     1.016095   0.412125   2.466  0.01368 *  
## sikap_num              0.968861   0.247448   3.915 9.02e-05 ***
## pendapatan_tahunan     0.014612   0.004499   3.248  0.00116 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 260.19  on 199  degrees of freedom
## Residual deviance: 225.13  on 195  degrees of freedom
## AIC: 235.13
## 
## Number of Fisher Scoring iterations: 4

Baseline yang Dipakai

Baseline adalah kategori referensi untuk variabel kategorikal (yang tidak muncul dalam output koefisien) atau nilai referensi untuk variabel numerik (biasanya 0).
- Variabel lokasi_hunian: Baseline adalah “rural”, karena koefisien hanya ditampilkan untuk “suburban” dan “urban”.
- Variabel sikap_num: Karena ini variabel numerik (ordinal yang dikonversi), tidak ada baseline kategorikal, tetapi kita bisa menganggap sikap_num = 0 sebagai titik referensi.
- Variabel pendapatan_tahunan: Juga numerik, sehingga referensinya adalah perubahan per unit dari nilai 0 atau nilai saat ini.

Interpretasi Setiap Koefisien

Berikut adalah penjelasan untuk setiap koefisien dalam model, termasuk log odds, p-value, dan odds ratio:

(Intercept): -5.308233

  • Log Odds: Log odds untuk membeli ketika semua prediktor berada pada baseline (lokasi_hunian = rural, sikap_num = 0, pendapatan_tahunan = 0) adalah -5.308233. Ini menunjukkan peluang membeli sangat rendah pada kondisi baseline.
  • Odds Ratio: exp(-5.308233) ≈ 0.0049, artinya peluang membeli hanya 0.49% pada baseline.
  • P-value: 7.14e-08 (sangat signifikan, p < 0.001), menunjukkan intercept berbeda secara signifikan dari nol.

lokasi_huniansuburban: 0.292033

  • Log Odds: Dibandingkan dengan “rural”, log odds untuk membeli di “suburban” meningkat sebesar 0.292033.
  • Odds Ratio: exp(0.292033) ≈ 1.339, artinya peluang membeli 1.339 kali lebih tinggi di “suburban” dibandingkan “rural”.
  • P-value: 0.48525 (tidak signifikan, p > 0.05), menunjukkan perbedaan ini tidak signifikan secara statistik.

lokasi_hunianurban: 1.016095

  • Log Odds: Dibandingkan dengan “rural”, log odds untuk membeli di “urban” meningkat sebesar 1.016095.
  • Odds Ratio: exp(1.016095) ≈ 2.762, artinya peluang membeli 2.76 kali lebih tinggi di “urban” dibandingkan “rural”.
  • P-value: 0.01368 (signifikan, p < 0.05), menunjukkan efek signifikan dari lokasi “urban”.

sikap_num: 0.968861

  • Log Odds: Setiap peningkatan satu unit pada sikap_num meningkatkan log odds untuk membeli sebesar 0.968861.
  • Odds Ratio: exp(0.968861) ≈ 2.635, artinya setiap peningkatan satu unit pada sikap_num meningkatkan peluang membeli 2.635 kali.
  • P-value: 9.02e-05 (sangat signifikan, p < 0.001), menunjukkan sikap_num adalah prediktor kuat.

pendapatan_tahunan: 0.014612

  • Log Odds: Setiap peningkatan satu unit pada pendapatan_tahunan (juta rupiah) meningkatkan log odds untuk membeli sebesar 0.014612.
  • Odds Ratio: exp(0.014612) ≈ 1.015, artinya setiap peningkatan satu unit pada pendapatan tahunan (juta rupiah) meningkatkan peluang membeli sebesar 1.5%.
  • P-value: 0.00116 (sangat signifikan, p < 0.01), menunjukkan efek signifikan.

Signifikansi Model

Variabel yang signifikan (p < 0.05) dalam meningkatkan peluang membeli adalah:
- lokasi_hunianurban (p = 0.01368)
- sikap_num (p = 9.02e-05)
- pendapatan_tahunan (p = 0.00116)

Variabel lokasi_huniansuburban (p = 0.48525) tidak signifikan, artinya tidak ada bukti kuat bahwa lokasi “suburban” berbeda dari “rural” dalam memengaruhi keputusan membeli.

Interpretasi Goodness-of-Fit

  • Null Deviance: 260.19 (deviance model tanpa prediktor, hanya intercept) dengan 199 derajat kebebasan.
  • Residual Deviance: 225.13 (deviance model saat ini) dengan 195 derajat kebebasan.
    • Penurunan deviance (260.19 - 225.13 = 35.06) menunjukkan model dengan prediktor lebih baik daripada model null.
  • AIC: 235.13 (Akaike Information Criterion), semakin rendah nilainya, semakin baik modelnya (dalam konteks ini, model numerik mempunyai AIC yang lebih kecil dibandingkan model nominal).

Analisis: Model ini lebih baik daripada model tanpa prediktor, tetapi residual deviance (225.13) yang masih besar menunjukkan ada variasi dalam data yang belum dijelaskan. Model mungkin bisa ditingkatkan dengan menambah prediktor atau interaksi antar variabel.

Kesimpulan Praktis

  • Faktor yang Meningkatkan Peluang Membeli:
    • Lokasi “urban”: Individu di area urban memiliki peluang membeli 2.76 kali lebih tinggi dibandingkan rural.
    • Tingkat kepedulian (sikap_num): Semakin tinggi nilai sikap_num, peluang membeli meningkat signifikan (2.635 kali per unit).
    • Pendapatan tahunan: Pendapatan yang lebih tinggi meningkatkan peluang membeli (1.5% per unit).
  • Faktor Tidak Signifikan: Lokasi “suburban” tidak berbeda signifikan dari “rural”.

Goodness-of-fit model

r2_mcfadden_model1 <- 1 - (model_nominal$deviance / model_nominal$null.deviance)
r2_mcfadden_model2 <- 1 - (model_num$deviance / model_num$null.deviance)
list(
McFadden_R2_Nominal = r2_mcfadden_model1,
McFadden_R2_Numerik = r2_mcfadden_model2
)
## $McFadden_R2_Nominal
## [1] 0.1353641
## 
## $McFadden_R2_Numerik
## [1] 0.1347507

Nilai R-squared McFadden digunakan untuk mengukur Goodness-of-fit, semakin besar nilainya maka semakin baik model dalam memprediksi.

Visualisasi Prediksi

Plot Model Nominal

library(ggplot2)
data_nominal$probabilitas <- predict(model_nominal, type = "response")
ggplot(data_nominal, aes(x = pendapatan_tahunan, y = probabilitas, color = sikap)) +
  geom_point(alpha = 0.6) +  
  labs(
    title = "Probabilitas Prediksi Membeli Kendaraan Listrik",
    subtitle = "Berdasarkan Pendapatan Tahunan, Lokasi Hunian, dan Sikap",
    x = "Pendapatan Tahunan (juta Rp)",
    y = "Probabilitas Membeli",
    color = "Sikap"
  ) +
  theme_minimal() +  
  theme(legend.position = "bottom") 

Plot Model Numerik

library(ggplot2)
data_numeric$probabilitas <- predict(model_num, type = "response")
ggplot(data_numeric, aes(x = pendapatan_tahunan, y = probabilitas, color = sikap_num)) +
  geom_point(alpha = 0.6) + 
  labs(
    title = "Probabilitas Prediksi Membeli Kendaraan Listrik",
    subtitle = "Berdasarkan Pendapatan Tahunan, Lokasi Hunian, dan Sikap Numerik",
    x = "Pendapatan Tahunan (juta Rp)",
    y = "Probabilitas Membeli",
    color = "Sikap Numerik"
  ) +
  scale_color_gradient(low = "blue", high = "red") +  
  theme_minimal() + 
  theme(legend.position = "bottom")  

10.4 Kesimpulan

  • Kedua model menunjukkan bahwa lokasi urban, sikap peduli (khususnya sangat peduli), dan pendapatan tahunan meningkatkan peluang membeli.

  • Model Numeric lebih sederhana (1 koefisien untuk sikap_num vs. 2 koefisien untuk sikapcukup_peduli dan sikapsangat_peduli), sehingga lebih mudah diinterpretasikan jika hubungan antar tingkat sikap dianggap linier.

  • Model Nominal lebih fleksibel karena tidak mengasumsikan linearitas, tetapi lebih kompleks dengan lebih banyak parameter.

11 Pemilihan Model Regresi Logistik dan Evaluasi

11.1 Pendekatan

Dalam analisis data, terdapat dua pendekatan utama untuk membangun model statistik: pendekatan konfirmatori dan eksploratori. Pendekatan konfirmatori bertujuan menguji hipotesis atau teori yang telah ditetapkan sebelumnya dengan menggunakan data untuk memvalidasi model tertentu. Sebaliknya, pendekatan eksploratori berfokus pada penemuan pola atau hubungan baru dalam data tanpa asumsi awal yang ketat. Kedua pendekatan ini memiliki teknik yang berbeda dan saling melengkapi, memungkinkan analis untuk memahami data secara mendalam dan membangun model yang sesuai dengan tujuan penelitian.

11.1.1 Confirmatory

Pendekatan Confirmatory dalam Regresi Logistik

Pendekatan Confirmatory adalah strategi pemilihan model yang bertujuan untuk menguji hipotesis spesifik atau teori tertentu dengan memilih prediktor berdasarkan pengetahuan domain, teori, atau penelitian sebelumnya, bukan hanya berdasarkan data itu sendiri. Berbeda dengan pendekatan eksplorasi (seperti forward/backward selection), pendekatan ini lebih terarah dan berfokus pada konfirmasi hubungan yang telah dihipotesiskan.

Karakteristik Pendekatan Confirmatory

  • Prediktor dipilih berdasarkan teori atau bukti sebelumnya, bukan hanya karena signifikansi statistik.

  • Ada hipotesis yang jelas tentang hubungan antara prediktor dan variabel dependen.

  • Penekanan pada interpretasi koefisien untuk mengkonfirmasi atau menolak hipotesis, bukan hanya memaksimalkan kecocokan model.

Contoh penggunaan :

Misalnya, teori menyatakan bahwa lokasi tempat tinggal dan sikap terhadap lingkungan memengaruhi probabilitas seseorang membeli kendaraan listrik. Maka model logistik dibangun langsung dengan lokasi_hunian dan sikap sebagai prediktor, lalu diuji untuk memverifikasi apakah sikap memberikan kontribusi signifikan.

11.1.2 Exploratory

Pendekatan Exploratory dalam Regresi Logistik

Pendekatan Explanatory adalah strategi pemilihan model yang bertujuan untuk menjelaskan hubungan antara prediktor dan variabel dependen berdasarkan data yang tersedia, dengan fokus pada eksplorasi pola dan hubungan yang mungkin belum dihipotesiskan sebelumnya. Pendekatan ini lebih fleksibel dan berorientasi pada penemuan hubungan baru dalam data.

Karakteristik Pendekatan Exploratory:

  • Prediktor dipilih berdasarkan data, sering kali melalui uji statistik seperti signifikansi atau kontribusi terhadap kecocokan model.

  • Tidak ada hipotesis spesifik yang ketat; tujuannya adalah menemukan prediktor yang paling menjelaskan variasi dalam data.

  • Penekanan pada maksimalisasi kecocokan model (misalnya, meminimalkan deviance atau AIC) untuk mengidentifikasi faktor-faktor penting.

Contoh penggunaan: Misalnya, untuk memahami faktor yang memengaruhi pembelian kendaraan listrik, model logistik awal dibuat dengan memasukkan berbagai prediktor seperti lokasi_hunian, sikap, pendapatan_tahunan, dan usia. Kemudian, prediktor yang tidak signifikan (misalnya, usia) dihapus melalui proses seperti backward elimination untuk menemukan model yang paling menjelaskan data.

Pembeda Tujuan Pendekatan Confirmatory dan Exploratory:

  • Confirmatory: Bertujuan untuk menguji hipotesis spesifik berdasarkan teori atau penelitian sebelumnya, memverifikasi apakah prediktor tertentu (misalnya, lokasi_hunian dan sikap) benar-benar memengaruhi hasil (contoh: pembelian kendaraan listrik) sesuai hipotesis.

  • Exploratory: Bertujuan untuk menemukan hubungan baru atau pola dalam data tanpa hipotesis awal yang ketat, dengan memilih prediktor (misalnya, lokasi_hunian, sikap, pendapatan_tahunan, usia) yang paling menjelaskan variasi data melalui metode seperti forward/backward selection.

11.2 Metode Stepwise : Forward, Backward, dan Kedua Arah

Metode forward selection, backward elimination, dan stepwise selection adalah pendekatan eksplorasi untuk pemilihan model regresi logistik yang bertujuan memilih kombinasi prediktor terbaik berdasarkan data. Forward selection dimulai dari model kosong dan menambahkan prediktor satu per satu yang paling meningkatkan kecocokan model (misalnya, berdasarkan AIC atau p-value). Backward elimination memulai dari model penuh dengan semua prediktor dan menghapus prediktor yang paling tidak signifikan secara bertahap. Stepwise selection menggabungkan keduanya, menambahkan atau menghapus prediktor berdasarkan kriteria kecocokan, seperti AIC, untuk menemukan model optimal.

  • Forward Selection: Mulai dari model kosong, tambahkan prediktor satu per satu yang paling meningkatkan kecocokan model.
  • Backward Elimination: Mulai dari model penuh, hapus prediktor satu per satu yang paling tidak signifikan.
  • Stepwise Selection: Gabungan forward dan backward, tambah atau hapus prediktor berdasarkan kriteria kecocokan model.

Simulasi :

library(knitr)
## Warning: package 'knitr' was built under R version 4.3.3
library(dplyr)
library(ggplot2)
library(MASS)
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
library(caret)
## Warning: package 'caret' was built under R version 4.3.3
## Loading required package: lattice
## Warning: package 'lattice' was built under R version 4.3.3
library(pROC)
## Warning: package 'pROC' was built under R version 4.3.3
## Type 'citation("pROC")' for a citation.
## 
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
## 
##     cov, smooth, var
library(DescTools)
## Warning: package 'DescTools' was built under R version 4.3.3
## 
## Attaching package: 'DescTools'
## The following objects are masked from 'package:caret':
## 
##     MAE, RMSE
set.seed(456)
n <- 250
var1 <- rnorm(n, mean = 0, sd = 1)
var2 <- rbinom(n, size = 1, prob = 0.4)
var3 <- rnorm(n, mean = 1, sd = 2)
logit <- -1 + 0.8 * var1 - 1.5 * var2 + 0.3 * var3
prob <- 1 / (1 + exp(-logit))
response <- rbinom(n, 1, prob)
data <- data.frame(outcome = as.factor(response), var1, var2, var3)
head(data)

Pembentukan Model :

model <- glm(outcome ~ var1 + var2 + var3, data = data, family = binomial)
summary(model)
## 
## Call:
## glm(formula = outcome ~ var1 + var2 + var3, family = binomial, 
##     data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -1.19331    0.25785  -4.628 3.69e-06 ***
## var1         0.89730    0.18138   4.947 7.53e-07 ***
## var2        -1.02715    0.33040  -3.109  0.00188 ** 
## var3         0.35086    0.08702   4.032 5.53e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 296.48  on 249  degrees of freedom
## Residual deviance: 239.45  on 246  degrees of freedom
## AIC: 247.45
## 
## Number of Fisher Scoring iterations: 5

Pemilihan Metode (Forward, Backward, Kedua Arah)

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

Dapat dilihat masing - masing nilai AIC dari ketiga metode (Forward, Backward, Kedua Arah) adalah sama.

11.3 Evaluasi Model : Kurva ROC dan AUC

ROC (Receiver Operating Characteristic)

ROC (Receiver Operating Characteristic) adalah plot visual dan analitis yang digunakan untuk mengevaluasi kinerja model klasifikasi, seperti regresi logistik, dengan memplot True Positive Rate (TPR) dengan False Positive Rate (FPR) pada berbagai nilai ambang batas (cut-off) probabilitas prediksi. ROC membantu menilai seberapa baik model membedakan antara dua kelas (misalnya, membeli vs tidak membeli kendaraan listrik).

Komponen Grafik ROC

  1. True Positive Rate (TPR) / Sensitivitas:
    • Definisi: Proporsi kasus positif (misalnya, membeli = 1) yang benar-benar diprediksi positif oleh model.
    • Rumus: TPR = TP / (TP + FN), di mana:
      • TP (True Positive): Jumlah kasus positif yang diprediksi benar.
      • FN (False Negative): Jumlah kasus positif yang diprediksi salah sebagai negatif.
    • Sumbu Y: TPR diplot pada sumbu vertikal, dengan rentang 0 hingga 1.
  2. False Positive Rate (FPR):
    • Definisi: Proporsi kasus negatif (misalnya, tidak membeli = 0) yang salah diprediksi sebagai positif.
    • Rumus: FPR = FP / (FP + TN), di mana:
      • FP (False Positive): Jumlah kasus negatif yang diprediksi salah sebagai positif.
      • TN (True Negative): Jumlah kasus negatif yang diprediksi benar.
    • Sumbu X: FPR diplot pada sumbu horizontal, dengan rentang 0 hingga 1.
  3. Titik pada Kurva:
    • Setiap titik pada kurva ROC mewakili pasangan (FPR, TPR) untuk ambang batas probabilitas tertentu (cut-off, misalnya, 0.5 atau 0.3).
    • Ambang batas adalah nilai probabilitas di mana model mengklasifikasikan prediksi sebagai positif (1) atau negatif (0).

Pergerakan Kurva (Cut-off)

  • Cut-off: Nilai ambang batas probabilitas yang menentukan apakah prediksi diklasifikasikan sebagai positif atau negatif. Misalnya, jika cut-off = 0.5, prediksi dengan probabilitas ≥ 0.5 diklasifikasikan sebagai 1.
  • Efek Perubahan Cut-off:
    • Cut-off rendah (misalnya, 0.1): Lebih banyak observasi diklasifikasikan sebagai positif, meningkatkan TPR tetapi juga FPR (kurva bergerak ke kanan atas).
    • Cut-off tinggi (misalnya, 0.9): Lebih sedikit observasi diklasifikasikan sebagai positif, menurunkan TPR dan FPR (kurva bergerak ke kiri bawah).
    • Pergerakan Kurva: Saat cut-off bervariasi dari 0 ke 1, kurva ROC menelusuri semua kombinasi TPR dan FPR, membentuk garis dari (0,0) ke (1,1).

Kurva ROC Ideal

  • Kurva ideal mencapai TPR = 1 dan FPR = 0 pada cut-off tertentu, membentuk sudut tajam di kiri atas grafik (titik (0,1)).
  • Ini berarti model mampu memisahkan kelas positif dan negatif dengan sempurna (tidak ada kesalahan klasifikasi).
  • Secara grafis, kurva mendekati sudut kiri atas dan memiliki AUC (Area Under the Curve) = 1.
  • Kurva Acak:
    • Jika model tidak lebih baik dari tebakan acak, kurva ROC akan mendekati garis diagonal (dari (0,0) ke (1,1)), dengan AUC = 0.5.
    • Ini menunjukkan model tidak memiliki kemampuan diskriminasi.

Contoh dalam Konteks

Dalam kasus adopsi kendaraan listrik:

  • Model regresi logistik (misalnya, model_num) menghasilkan probabilitas prediksi (probabilitas).

  • ROC diplot dengan menghitung TPR dan FPR untuk berbagai cut-off (misalnya, probabilitas 0.3, 0.5, 0.7).

  • Jika kurva ROC mendekati sudut kiri atas dan AUC mendekati 1, maka model sangat baik dalam membedakan antara individu yang membeli dan tidak membeli kendaraan listrik.

Kesimpulan

ROC adalah alat visual untuk menilai performa klasifikasi model dengan memplot TPR vs FPR pada berbagai cut-off. Kurva ideal mencapai TPR = 1 dan FPR = 0 (AUC = 1), sementara pergerakan kurva dipengaruhi oleh perubahan cut-off, yang memengaruhi trade-off antara sensitivitas dan spesifisitas.

Contoh (Dari Data Simulasi)

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

Memilih Threshold yang Baik

AUC (Area Under Curve)

Mengukur luas di bawah kurva ROC, menunjukkan kemampuan model membedakan kelas positif dan negatif. AUC mengukur luas di bawah kurva ROC, dengan rentang 0 hingga 1:

  • AUC = 1: Model sempurna, membedakan kelas positif dan negatif tanpa kesalahan.
  • AUC = 0.5: Model tidak lebih baik dari tebakan acak.

AUC > 0.5: Model memiliki kemampuan diskriminasi yang baik, semakin mendekati 1 semakin baik.

Perhitungan (Dari Data Simulasi)

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

Terlihat dari nilai AUC-nya bahwa model sudah lebih baik dari tebakan acak walau belum mendekati sempurna.

11.4 Evaluasi Model : Pseudo R-squared

Pseudo R-Squared adalah ukuran yang digunakan untuk mengevaluasi kecocokan model regresi logistik, mirip dengan R² pada regresi linier, tetapi disesuaikan untuk model dengan variabel dependen biner. Berbeda dengan R² linier, pseudo R-squared tidak mengukur proporsi variansi yang dijelaskan secara langsung, melainkan memberikan indikasi seberapa baik model menjelaskan data dibandingkan model dasar (null model).

McFadden R² fokus pada pengurangan deviance, memberikan ukuran proporsi “penjelasan” model dibandingkan model tanpa prediktor. Nilai yang lebih tinggi (R² >0.2) menunjukkan model yang cukup baik, terutama dalam konteks regresi logistik di mana data biner sering sulit dijelaskan sepenuhnya.

Cox & Snell R² menggunakan rasio likelihood untuk mengevaluasi peningkatan kecocokan, tetapi cenderung lebih kecil dan kurang interpretable dibandingkan McFadden karena batasan matematisnya.

Perhitungan (Dari Data Simulasi)

PseudoR2(step_both, which = c("CoxSnell", "McFadden"))
##  CoxSnell  McFadden 
## 0.2039576 0.1923447

11.5 Evaluasi Model : Sensitifitas dan Spesifisitas

Sensitifitas dan Spesifisitas

  • Sensitivity (True Positive Rate): Proporsi kasus positif (misalnya, membeli kendaraan listrik) yang benar diprediksi positif oleh model. Rumus: TP / (TP + FN).
  • Specificity (True Negative Rate): Proporsi kasus negatif (misalnya, tidak membeli) yang benar diprediksi negatif. Rumus: TN / (TN + FP).

Contoh (Dari Data Simulasi)

pred_class <- ifelse(pred_prob >= 0.5, 1, 0)
conf_matrix <- confusionMatrix(factor(pred_class), data$outcome, positive = "1")
conf_matrix$byClass[c("Sensitivity", "Specificity")]
## Sensitivity Specificity 
##   0.4000000   0.9277778

Dari semua positif aktual, hanya 40% yang berhasil diprediksi positif oleh model dan banyak false negative. Dari semua negatif aktual, 92.78% berhasil diprediksi negatif oleh model sehingga false positive tergolong sedikit.

Kesimpulan (Dari Data Simulasi)

  • Model sangat baik dalam mendeteksi yang negatif (spesifik tinggi)

  • Tapi kurang bagus dalam mendeteksi yang positif (sensitivitas rendah).

11.6 Evaluasi Model : AIC dan Deviance

AIC (Akaike Information Criterion) mengukur kualitas model dengan menyeimbangkan kecocokan data (deviance) dan kompleksitas (jumlah parameter), di mana nilai lebih rendah menunjukkan model lebih baik. Deviance mengukur seberapa baik model menjelaskan data dibandingkan model jenuh (sempurna), dengan nilai lebih rendah menunjukkan kecocokan yang lebih baik. Idealnya, AIC minimal untuk menghindari overfitting, sementara deviance rendah menunjukkan model mendekati data aktual, tetapi harus diimbangi dengan jumlah prediktor agar tidak terlalu kompleks.

Hubungan keduanya adalah bahwa penurunan deviance (meningkatkan kecocokan) dengan penambahan prediktor dapat meningkatkan AIC jika kompleksitas model terlalu tinggi, sehingga memilih model terbaik membutuhkan trade-off. Dalam praktiknya, pilih model dengan AIC terendah yang juga menunjukkan penurunan deviance signifikan untuk memastikan efisiensi dan keandalan prediksi.

Contoh (Dari Data Simulasi)

Pembentukan Model

model1 <- glm(outcome ~ var1, data = data, family = binomial)
model2 <- glm(outcome ~ var1 + var2, data = data, family = binomial)
model3 <- glm(outcome ~ var1 + var2 + var3, data = data, family = binomial)

Perbandingan AIC dan Deviance

perbandingan <- 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))
)
perbandingan

Kesimpulan

  • Model Terbaik: Model 2, karena memiliki AIC terendah (264.2116), menunjukkan keseimbangan terbaik antara kecocokan dan kompleksitas.

  • Deviance: Penurunan dari Model 1 ke Model 3 menunjukkan prediktor tambahan meningkatkan kecocokan, tetapi Model 3 terlalu kompleks (AIC naik).

  • Gunakan Model 2 untuk prediksi atau analisis, karena lebih efisien dan cukup menjelaskan data tanpa overfitting.

11.7 Prinsip Parsimony

Prinsip Parsimony adalah konsep yang menyatakan bahwa dari beberapa model yang menjelaskan data dengan baik, model yang paling sederhana (dengan jumlah parameter atau asumsi paling sedikit) harus dipilih, kecuali ada bukti kuat mendukung model yang lebih kompleks.

11.8 Evaluasi Model : Likelihood-Ratio Test

Likelihood-Ratio Test (LRT) dalam konteks regresi logistik digunakan untuk membandingkan dua model dengan tujuan untuk menguji apakah penambahan prediktor secara signifikan meningkatkan kecocokan model.:

  • Model penuh (full model): dengan semua prediktor.

  • Model terbatas (reduced model): dengan sebagian prediktor (biasanya hanya konstanta).

Rumus Likelihood Ratio\[ \text{LR} = -2 \left( \log L_0 - \log L_1 \right) \]

  • \(\log L_0\): log-likelihood model terbatas

  • \(\log L_1\): log-likelihood model penuh

  • LR mengikuti distribusi \(\chi^2\) (chi-square) dengan derajat bebas = selisih jumlah parameter antara dua model.

Contoh Perhitungan (Dari Data Simulasi)

model <- glm(outcome ~ var1 + var2 + var3, data = data, family = binomial)
model_null <- glm(outcome ~ 1, data = data, family = binomial)
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

11.9 Evaluasi Model : Precision-Recall Curve (PR Curve)

Recall (Sensitivity / TPR) \[ \text{Recall} = \frac{\text{TP}}{\text{TP + FN}} \] Artinya: Dari semua positif aktual, berapa persen yang berhasil diprediksi positif.

Precision (PPV) \[ \text{Precision} = \frac{\text{TP}}{\text{TP + FP}} \] Artinya: Dari semua yang diprediksi positif, berapa persen yang benar-benar positif.

Kurva PR

  • Grafik: Precision (y-axis) vs Recall (x-axis).

  • Dibuat dengan menghitung precision dan recall di berbagai threshold probabilitas.

Area Under PR Curve (AUPRC)

  • Luas area di bawah PR curve.

  • Semakin besar AUPRC (maks = 1), semakin baik model mendeteksi positif dengan minim false positive.

  • Nilai baseline AUPRC = proporsi positif di data.

Perbandingan dengan ROC Curve

Aspek PR Curve ROC Curve
Fokus Kinerja pada kelas positif Keseimbangan antara positif & negatif
Y-Axis Precision True Positive Rate (Recall)
X-Axis Recall False Positive Rate (FPR)
Cocok untuk Data imbalanced (positif jarang) Data seimbang
Baseline Proporsi positif 0.5

Visualisasi PR Curve

library(PRROC)
## Warning: package 'PRROC' was built under R version 4.3.3
## Loading required package: rlang
set.seed(123)
x1 <- rnorm(200)
x2 <- rbinom(200, 1, 0.5)
x3 <- rnorm(200)
lin_pred <- -1 + 1.5 * x1 - 0.7 * x2 + 0.6 * x3
p <- 1 / (1 + exp(-lin_pred))
y <- rbinom(200, 1, p)
data <- data.frame(y = y, x1, x2, x3)
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)

Precision-Recall Curve (PR Curve) cocok diimplementasikan saat:

  • Data imbalanced (kelas positif jauh lebih sedikit dari negatif)

  • Fokus pada kinerja deteksi kelas positif

  • Ingin meminimalkan false positive dan false negative

  • ROC Curve tampak bagus, tapi model tetap gagal menangkap positif dengan baik

ROC Curve cocok digunakan pada beberpa aplikasi nyata, contoh : deteksi penipuan , spam filter.

12 Regresi Logistik Multinomial

12.1 Distribusi Multinomial

Distribusi Multinomial adalah generalisasi distribusi binomial untuk lebih dari dua kategori hasil. Ini menggambarkan probabilitas kejadian dalam \(n\) percobaan independen, di mana setiap percobaan memiliki \(k\) kategori hasil dengan probabilitas masing-masing \(p_1, p_2, \dots, p_k\), dan \(\sum p_i = 1\).

Probabilitas

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

Keterangan:
\(n\): Jumlah total percobaan.
\(x_i\): Jumlah kejadian untuk kategori \(i\), dengan \(\sum x_i = n\).
\(p_i\): Probabilitas kategori \(i\).

Contoh Kasus Sederhana
Sebuah dadu dengan tiga sisi (merah, biru, hijau) dilempar 5 kali. Probabilitas munculnya masing-masing sisi adalah \(p_{\text{merah}} = 0.4\), \(p_{\text{biru}} = 0.3\), \(p_{\text{hijau}} = 0.2\). Hitung probabilitas mendapatkan 2 merah, 2 biru, dan 1 hijau.

Perhitungan R

n <- 5
x <- c(2, 2, 1)  # merah, biru, hijau
p <- c(0.4, 0.3, 0.2)
prob <- dmultinom(x, size = n, prob = p)
print(prob)
## [1] 0.1463192

Jadi, probabilitas mendapatkan sisi dadu dengan komposisi : 2 merah, 2 biru, dan 2 hijau adalah 0.14631 atau 14.63%.

12.2 Regresi Logistik Multinomial

Regresi Logistik Multinomial adalah model statistik untuk memprediksi probabilitas hasil kategori dengan lebih dari dua kelas (misalnya, A, B, C) berdasarkan variabel independen. Model ini memperluas regresi logistik biner dengan menggunakan fungsi eksponensial ternormalisasi untuk menghitung probabilitas setiap kelas relatif terhadap kelas referensi.

Perbedaan dengan Regresi Logistik Biasa:

  • Jumlah Kategori: Regresi logistik biasa (biner) hanya menangani dua kategori hasil (misalnya, ya/tidak), sedangkan multinomial menangani lebih dari dua kategori.

  • Fungsi logit:

    \[ \log\left(\frac{P(y_i = k | x_1, x_2, \dots, x_p)}{P(y_i = K | x_1, x_2, \dots, x_p)}\right) = \beta_{k0} + \beta_{k1}x_1 + \beta_{k2}x_2 + \dots + \beta_{kp}x_p, \quad k = 1, 2, \dots, K-1 \]

    Keterangan:
    \(P(y_i = k | x_1, x_2, \dots, x_p)\): Probabilitas kelas \(k\) untuk observasi dengan variabel independen \(x_1, x_2, \dots, x_p\).
    \(P(y_i = K | x_1, x_2, \dots, x_p)\): Probabilitas kelas referensi \(K\).
    \(\beta_{k0}\): Intersep untuk kelas \(k\).
    \(\beta_{k1}, \beta_{k2}, \dots, \beta_{kp}\): Koefisien untuk variabel independen \(x_1, x_2, \dots, x_p\) pada kelas \(k\).
    \(K\): Jumlah kelas.

  • Interpretasi: Dalam regresi multinomial, koefisien diinterpretasikan relatif terhadap kelas referensi, bukan hanya antara dua kelas seperti pada regresi biner.

12.2.1 Baseline-category logit model

Model baseline-category logit merupakan bentuk regresi logistik yang digunakan untuk variabel respon kategorik dengan lebih dari dua kategori yang bersifat nominal. Dalam model ini, salah satu kategori dipilih sebagai kategori acuan (baseline), dan kategori lainnya dibandingkan terhadap acuan tersebut dalam bentuk fungsi logit.

Untuk setiap kategori \(j = 1, 2, \dots J-1\), model dituliskan sebagai:

\[\log\left( \frac{\pi_j(\mathbf{x})}{\pi_J(\mathbf{x})} \right) = \alpha_j + \boldsymbol{\beta}_j^\top \mathbf{x}\]

12.2.2 Estimasi Parameter

Dalam regresi logistik multinomial (terutama model baseline-category logit), parameter diestimasi menggunakan metode Maximum Likelihood Estimation (MLE).

Tujuannya adalah mencari nilai parameter (intersep dan koefisien regresi) yang memaksimalkan kemungkinan (likelihood) data yang diamati.

Estimasi Parameter dengan Maximum Likelihood

Parameter \(\beta_j\) (untuk j = 1, 2, …, J-1) diestimasi dengan memaksimalkan fungsi log-likelihood. Untuk n observasi, dengan \(y_i\) sebagai kategori untuk observasi ke-i dan \(x_i\) sebagai vektor kovariat, fungsi likelihood adalah:

\[L(\beta) = \prod_{i=1}^n \prod_{j=1}^J [P(Y_i = j | x_i)]^I(y_i = j)\]

12.3 Contoh Kasus

Regresi Logistik Multinomial dapat digunakan untuk memodelkan probabilitas memilih salah satu merek (Coca-Cola, Pepsi, atau Sprite) sebagai variabel dependen kategorikal dengan tiga kelas. Variabel independennya adalah:

  • Usia: Numerik (misalnya, dalam tahun).

  • Frekuensi konsumsi minuman per minggu: Numerik (misalnya, jumlah minuman ringan yang dikonsumsi per minggu).

  • Jenis kelamin: Kategorikal (misalnya, laki-laki/perempuan).

12.4 Simulasi Data

set.seed(456)
n <- 100  
usia <- sample(20:40, 2*n, replace = TRUE)
frekuensi <- sample(1:7, 2*n, replace = TRUE)
jenkel <- rep(c("Laki-laki", "Perempuan"), each = n)

prob_laki <- c(0.5, 0.3, 0.2)  # Coca-Cola, Pepsi, Sprite
prob_perempuan <- c(0.2, 0.4, 0.4)

# Simulasi merek
merek_laki <- sample(c("Coca-Cola", "Pepsi", "Sprite"), n, replace = TRUE, prob = prob_laki)
merek_perempuan <- sample(c("Coca-Cola", "Pepsi", "Sprite"), n, replace = TRUE, prob = prob_perempuan)
merek <- c(merek_laki, merek_perempuan)

# Gabungkan data
data_minuman <- data.frame(
  Usia = usia,
  Frekuensi = frekuensi,
  Jenis_Kelamin = jenkel,
  Merek = merek
)
head(data_minuman)

12.5 Estimasi Model

library(nnet)
## Warning: package 'nnet' was built under R version 4.3.3
model_multi <- multinom(Merek ~ Usia + Frekuensi + Jenis_Kelamin, 
                          data = data_minuman)
## # weights:  15 (8 variable)
## initial  value 219.722458 
## iter  10 value 211.247374
## final  value 211.241812 
## converged
summary(model_multi)
## Call:
## multinom(formula = Merek ~ Usia + Frekuensi + Jenis_Kelamin, 
##     data = data_minuman)
## 
## Coefficients:
##        (Intercept)        Usia  Frekuensi Jenis_KelaminPerempuan
## Pepsi   -2.5675377 0.067150907 0.02927628              0.6677877
## Sprite  -0.9009979 0.009936885 0.09105572              0.8718447
## 
## Std. Errors:
##        (Intercept)       Usia  Frekuensi Jenis_KelaminPerempuan
## Pepsi    1.0409313 0.03225753 0.09056037              0.3744939
## Sprite   0.9417589 0.03000216 0.08360399              0.3479144
## 
## Residual Deviance: 422.4836 
## AIC: 438.4836

12.6 Nilai P-value dan Interpretasi

z <- summary(model_multi)$coefficients / summary(model_multi)$standard.errors
pval <- 2 * (1 - pnorm(abs(z)))
round(pval, 4)
##        (Intercept)   Usia Frekuensi Jenis_KelaminPerempuan
## Pepsi       0.0136 0.0374    0.7465                 0.0746
## Sprite      0.3387 0.7405    0.2761                 0.0122

Interpretasi:

• Baseline untuk merek minuman adalah “Coca-cola”

• Variabel signifikan : hanya variabel usia (berpengaruh terhadap peluang memilih Pepsi vs Coca-cola) dan jenis kelamin (berpengaruh terhadap peluang memilih Sprite vs Coca-Cola) yang signifikan.

12.7 Prediksi dan Validasi

data_minuman$Predicted <- predict(model_multi)
table(Predicted = data_minuman$Predicted, Actual = data_minuman$Merek)
##            Actual
## Predicted   Coca-Cola Pepsi Sprite
##   Coca-Cola        37    21     24
##   Pepsi             8     7     11
##   Sprite           21    29     42

12.8 Kesimpulan

Model ini tidak cukup bagus dalam memprediksi preferensi merek. Akurasi rendah, dan performa buruk pada kelas Pepsi menunjukkan model mungkin memerlukan lebih banyak data, variabel independen yang lebih relevan, atau penyesuaian parameter (misalnya, regularisasi atau penanganan ketidakseimbangan kelas).

13 Regresi Logistik Ordinal

Regresi Logistik Ordinal adalah metode regresi yang menganalisis hubungan antara satu atau lebih variabel prediktor dengan variabel respons ordinal, yaitu variabel kategori dengan urutan misal: Tingkat keparahan penyakit (ringan, sedang, parah) dll.

13.1 Konsep Cumulative Logit Model

Untuk variabel respons ordinal Y dengan \(J\) kategori berurutan \[(1,2,…\dots, J)\] model logit kumulatif didefinisikan sebagai:

\(\log \left( \frac{P(Y \leq j)}{P(Y > j)} \right) = \theta_j - \mathbf{x}^\top \boldsymbol{\beta}, \quad j = 1, 2, \dots, J - 1\)

  • \(P(Y \leq j)\): Probabilitas kumulatif bahwa respons berada pada kategori ke-j atau kurang

  • \(θj\): Intersep spesifik kategori (threshold/cutpoint)

  • \(\mathbf{x}\): Vektor variabel prediktor

  • \(\boldsymbol{\beta}\): Koefisien regresi (sama untuk semua j)

Konsep:

  • Model ini digunakan saat variabel dependen bersifat ordinal (terurut), seperti tingkat keparahan penyakit : ringan, sedang, dan parah.

  • Cumulative logit model menghitung logit dari probabilitas kumulatif untuk setiap batas kategori.

  • Koefisien \(\boldsymbol{\beta}\) menggambarkan pengaruh prediktor terhadap peluang naik ke kategori lebih tinggi.

Asumsi :

  • Proportional Odds Assumption (Asumsi peluang proporsional):
    Koefisien \(\boldsymbol{\beta}\) tetap sama untuk semua batas kumulatif j.
    Ini berarti efek prediktor bersifat konsisten lintas tingkat kategorisasi.

  • Independensi Observasi: Setiap observasi independen satu sama lain.

  • Kategori bersifat ordinal: Harus ada urutan yang jelas antara kategori (meski jarak antar kategori tidak harus sama).

13.2 Interpretasi Koefisien

Koefisien \(\boldsymbol{\beta}\) dalam regresi logistik ordinal menginterpretasikan arah dan kekuatan pengaruh variabel prediktor terhadap peluang berada pada kategori yang lebih tinggi dari variabel respons ordinal.

  • Jika \(\beta > 0\):
    Variabel prediktor meningkatkan peluang individu berada pada kategori yang lebih tinggi.
    (Efek positif terhadap “kenaikan level”).

  • Jika \(\beta < 0\):
    Variabel prediktor menurunkan peluang berada pada kategori lebih tinggi.
    (Efek negatif terhadap “kenaikan level”).

  • Jika \(\beta = 0\):
    Tidak ada pengaruh prediktor terhadap respons ordinal.

Besar Kecilnya Koefisien

  • Koefisien besar (absolut) maka pengaruh kuat terhadap perubahan tingkat kategori.

  • Koefisien kecil (mendekati nol) maka pengaruh lemah terhadap perubahan tingkat kategori.

Interpretasi Odds Ratio

Kadang koefisien diinterpretasi lewat odds ratio:

\(\text{Odds Ratio}\) = \(e^{\beta}\)

  • Jika \(e^{\beta} > 1\): peluang naik kategori meningkat

  • Jika \(e^{\beta} < 1\): peluang naik kategori menurun

13.3 Contoh Data

Data ini berisi 100 observaasi hasil simulasi dari penggolongan keparahan hipertensi yang terurut (ringan, sedang, parah) berdasarkan prediktor usia.
*Data bersifat dummy dan tidak bisa dikaitkan dengan kasus riil.

library(MASS)
set.seed(123)
n <- 100
usia <- round(runif(n, 20, 70))  
lp <- 0.08 * usia + rnorm(n, 0, 1)

# Potong menjadi 3 kategori ordinal
keparahan <- cut(lp,
                 breaks = quantile(lp, probs = c(0, 1/3, 2/3, 1)),
                 labels = c("Ringan", "Sedang", "Parah"),
                 include.lowest = TRUE,
                 ordered_result = TRUE)

# Gabung ke data frame
data_sim <- data.frame(keparahan, usia)
table(data_sim$keparahan)
## 
## Ringan Sedang  Parah 
##     34     33     33

13.4 Estimasi Model Ordinal

# Model regresi logistik ordinal (proportional odds)
model <- polr(keparahan ~ usia, data = data_sim, Hess = TRUE)
summary(model)
## Call:
## polr(formula = keparahan ~ usia, data = data_sim, Hess = TRUE)
## 
## Coefficients:
##       Value Std. Error t value
## usia 0.1373    0.02045   6.715
## 
## Intercepts:
##               Value  Std. Error t value
## Ringan|Sedang 4.9793 0.8574     5.8074 
## Sedang|Parah  7.4014 1.0652     6.9481 
## 
## Residual Deviance: 152.6484 
## AIC: 158.6484

13.5 Nilai P-value

# Hitung p-value manual
ctable <- coef(summary(model))
pval <- pnorm(abs(ctable[, "t value"]), lower.tail = FALSE) * 2
ctable <- cbind(ctable, "p value" = pval)
ctable
##                   Value Std. Error  t value      p value
## usia          0.1373114 0.02044961 6.714620 1.885565e-11
## Ringan|Sedang 4.9792755 0.85740445 5.807382 6.345711e-09
## Sedang|Parah  7.4013777 1.06523828 6.948096 3.702493e-12

Dari p-value nya, dapat disimpulkan bahwa usia secara signifikan meningkatkan kemungkinan seseorang mengalami keparahan hipertensi yang lebih tinggi. Model memperlihatkan bahwa dengan bertambahnya usia, risiko berpindah dari kategori “Ringan” ke “Sedang” atau dari “Sedang” ke “Parah” juga meningkat secara signifikan.

13.6 Prediksi Probabilitas

# Prediksi probabilitas untuk usia tertentu
newdata <- data.frame(usia = c(25, 40, 60))
pred <- predict(model, newdata = newdata, type = "probs")
hasil_prediksi <- data.frame(usia = newdata$usia, pred)
print(hasil_prediksi)
##   usia     Ringan    Sedang      Parah
## 1   25 0.82440644 0.1570442 0.01854941
## 2   40 0.37444867 0.4964495 0.12910183
## 3   60 0.03699037 0.2651126 0.69789699

Tabel di atas menunjukkan probabilitas masuknya individu dengan usia 25, 40, dan 60 ke dalam kategori keparahan penyakit hipertensi.

Terlihat dari tabel bahwa semakin muda usianya, semakin kecil peluangnya mendapatkan keparahan hipertensi pada tingkat “parah”.

13.7 Goodness-of-fit dan Proportional Odds

Model cumulative logit mengasumsikan bahwa pengaruh prediktor tetap konsisten di setiap batas kategori. Jika asumsi ini tidak terpenuhi, maka model alternatif seperti generalized ordinal model yang tidak mengharuskan proporsionalitas odds dapat dipertimbangkan.

13.8 Alternatif Model Ordinal

Alternatif dari regresi logistik ordinal (terutama proportional odds model) dapat digunakan ketika asumsi-asumsinya seperti proportional odds.

  1. Generalized Ordinal Logistic Regression
  • Cocok saat: Asumsi proportional odds tidak terpenuhi.
  • Keunggulan: Mengizinkan efek prediktor berbeda di setiap cutpoint.
  1. Partial Proportional Odds Model (PPOM)
  • Cocok saat: Hanya sebagian prediktor melanggar asumsi proportional odds.
  • Keunggulan: Fleksibel karena hanya prediktor tertentu yang diberi efek berbeda di tiap kategori.
  1. Adjacent Category Logistic Model
  • Cocok saat: Ingin memodelkan logit antara pasangan kategori berurutan.

  • Keunggulan: Memperhitungkan struktur ordinal secara eksplisit.

13.10 Asumsi Paralelisme

Asumsi paralelisme (proportional odds) menyatakan bahwa: pengaruh setiap prediktor terhadap peluang kumulatif dari kategori respon bersifat konstan (paralel) di seluruh batas antar kategori. Dengan kata lain, koefisien dari prediktor tidak berubah untuk setiap batas antar kategori respon ordinal.

Jika asumsi ini tidak terpenuhi, maka koefisien prediktor berubah-ubah tergantung batas kategori mana yang dibandingkan dan model ordinal logistic jadi tidak cocok sehingga dibutuhkan alternatif yaitu:

  • Generalized ordered logit model atau

  • Multinomial logit model (tanpa asumsi paralel)

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

Karena p > 0.05, maka tidak ada pelanggaran signifikan terhadap asumsi proportional odds.

14 Log Linear Model

Log-linear model merupakan salah satu pendekatan statistik yang sangat penting untuk menganalisis data kategorik, khususnya ketika data tersebut disajikan dalam bentuk tabel kontingensi multidimensi. Pendekatan digunakan untuk memahami pola hubungan dan interaksi antara dua atau lebih variabel kategorik yang tercermin dalam frekuensi atau jumlah kasus pada setiap kombinasi kategori. Berbeda dengan analisis tabel kontingensi biasa yang hanya mampu menguji asosiasi sederhana antar variabel (misalnya dengan uji chi-square), log-linear model menggunakan kerangka yang lebih fleksibel dan menyeluruh dengan memodelkan logaritma dari ekspektasi frekuensi sel sebagai fungsi linear dari parameter-parameter yang berkaitan dengan variabel-variabel dalam tabel tersebut. Dalam konteks ini, model log-linear dapat dianggap sebagai perluasan dari analisis tabel kontingensi, di mana bukan hanya asosiasi dua variabel yang dapat diuji, tetapi juga interaksi lebih dari dua variabel sekaligus, serta efek utama (main effect) masing-masing variabel.

Secara struktur, log-linear model membangun hubungan antara variabel-variabel kategorik dengan memodelkan logaritma dari ekspektasi frekuensi sel sebagai kombinasi linear dari parameter-parameter yang merepresentasikan efek utama dan interaksi antar variabel . Jika interaksi antar variabel tidak signifikan, maka parameter interaksi dapat dihilangkan, menandakan bahwa kedua variabel bersifat independen secara statistik. Namun, jika terdapat asosiasi atau ketergantungan, maka parameter interaksi menjadi sangat penting untuk dimasukkan ke dalam model. Untuk mencegah overparameterisasi dan memastikan identifikasi model, biasanya digunakan constraint seperti sum-to-zero, yaitu jumlah seluruh parameter efek utama sama dengan nol.

Tujuan utama dari log-linear model adalah untuk mengidentifikasi dan mengukur kekuatan asosiasi serta interaksi antara variabel-variabel kategorik dalam data, sekaligus menguji model mana yang paling sesuai dengan struktur data yang ada. Model ini sangat berguna dalam situasi di mana data berbentuk hitungan (count data), seperti jumlah kasus penyakit berdasarkan kelompok usia dan jenis kelamin, atau pola pembelian produk berdasarkan kategori demografi. Selain itu, log-linear model juga memungkinkan pengujian model-model hierarkis, mulai dari model dengan hanya efek utama hingga model dengan interaksi tingkat tinggi (saturated model), sehingga peneliti dapat memilih model yang paling parsimonious namun tetap mampu menjelaskan data dengan baik.

Dari sisi hubungan dengan Generalized Linear Model (GLM), log-linear model merupakan salah satu spesifikasi dari GLM dengan asumsi bahwa distribusi data mengikuti distribusi Poisson dan fungsi link yang digunakan adalah logaritma alami (natural log). Dengan kata lain, log-linear model memodelkan ekspektasi frekuensi sel sebagai fungsi eksponensial dari kombinasi linear parameter, sehingga sangat cocok untuk data hitungan yang tidak bernilai negatif dan sering kali bersifat skewed. Transformasi logaritma pada frekuensi sel juga membantu mengatasi masalah non-normalitas pada data count, sehingga hasil estimasi parameter mudah diinterpretasikan setelah dilakukan back-transformation.

Dalam praktiknya, log-linear model sangat bermanfaat di berbagai bidang seperti epidemiologi, ilmu sosial, pemasaran, dan biostatistika, di mana analisis hubungan dan interaksi antar kategori sangat krusial. Keunggulan utama model ini adalah kemampuannya untuk menangani data multi-dimensi dan memberikan pemahaman yang lebih mendalam tentang struktur asosiasi dalam data kategorik. Namun, penting untuk memperhatikan asumsi-asumsi dasar seperti distribusi Poisson dan kecukupan jumlah data pada setiap sel tabel, agar hasil analisis tetap valid dan dapat diandalkan.

Sebagai kesimpulan, log-linear model merupakan alat analisis statistik yang penting untuk memahami dan memodelkan hubungan kompleks antara variabel-variabel kategorik dalam data berbentuk tabel kontingensi. Dengan struktur model yang fleksibel, kemampuan menangani interaksi multi-variabel, serta keterkaitannya dengan kerangka GLM, log-linear model menjadi fondasi penting dalam analisis data kategorik, melebihi keterbatasan analisis tabel kontingensi konvensional dan membuka peluang eksplorasi pola-pola asosiasi yang lebih kompleks.

Perbedaan Tabel Kontingensi, Regresi Logistik, dan Model Log-linear

Aspek Tabel Kontingensi Regresi Logistik Model Log-Linear
Definisi Penyajian frekuensi gabungan dua atau lebih variabel kategorik dalam bentuk tabel Model statistik untuk memprediksi probabilitas kategori dependen berdasarkan variabel prediktor Model statistik untuk memodelkan log frekuensi sel tabel kontingensi sebagai fungsi efek variabel dan interaksinya
Tujuan Analisis Eksplorasi hubungan dan asosiasi antar variabel kategorik secara sederhana Mengestimasi hubungan antara variabel dependen kategorik dengan satu atau lebih variabel independen Menguji dan memodelkan asosiasi dan interaksi antar variabel kategorik dalam tabel multidimensi
Variabel Dependen Tidak ada variabel dependen khusus, hanya frekuensi sel Ada variabel dependen kategorik (biner atau multinomial) Tidak membedakan variabel dependen dan independen, semua variabel diperlakukan simetris
Jenis Data Data frekuensi atau count dalam tabel silang Data kategorik dengan variabel respon kategorik dan prediktor kategorik/kontinu Data count (frekuensi) dalam tabel kontingensi multidimensi
Asumsi Distribusi Tidak spesifik, hanya frekuensi observasi Variabel respon mengikuti distribusi Bernoulli atau multinomial Frekuensi sel mengikuti distribusi Poisson
Pendekatan Analisis Uji chi-square untuk independensi dan asosiasi Estimasi parameter dengan likelihood, uji signifikansi koefisien Estimasi parameter dengan likelihood, uji goodness-of-fit model
Interaksi Variabel Dilihat secara sederhana melalui tabel silang Dapat dimodelkan sebagai interaksi prediktor dalam model Interaksi antar variabel secara eksplisit dimodelkan sebagai parameter interaksi
Kelebihan Mudah dibuat dan dipahami, cocok untuk eksplorasi awal Memungkinkan prediksi probabilitas dan kontrol variabel lain Fleksibel untuk model multidimensi, menangani interaksi kompleks dan asosiasi

14.1 Tabel Kontingensi

Untuk mendalami model log linier, akan disimulasikan data pelatihan dan produktivitas pada karyawan di sebuah perusahaan. Data merupakan data tabel kontingensi 2x2 :

data_pelatihan <- array(
  c(75, 25, 40, 60), 
  dim = c(2, 2),
  dimnames = list(
    Pelatihan = c("Ikut", "Tidak Ikut"),
    Produktivitas = c("Tinggi", "Rendah")
  )
)
data_pelatihan
##             Produktivitas
## Pelatihan    Tinggi Rendah
##   Ikut           75     40
##   Tidak Ikut     25     60

14.2 Model Saturated

Model saturated mencakup semua efek utama (main effects) dan semua interaksi antar variabel sampai tingkat tertinggi.

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

Dengan:

\(\mu_{ij}\) adalah nilai harapan (mean) untuk sel pada kategori \(i\) dari faktor \(A\) dan \(j\) dari faktor \(B\).

\(\lambda\): konstanta (intersep)

\(\lambda^A_i\) : efek utama dari variabel A (i = 1,2)

\(\lambda^B_j\)​: efek utama: efek utama dari variabel B (j = 1,2)

\(\lambda^{AB}_{ij}\) : efek interaksi antara A dan B

Dari data tabel kontingensi yang sudah dibuat, maka didapatkan model saturatednya :

Contoh (Data Pelatihan)

library(MASS)
model_saturated <- loglm(~ Pelatihan * Produktivitas, data = data_pelatihan)
summary(model_saturated)
## Formula:
## ~Pelatihan * Produktivitas
## attr(,"variables")
## list(Pelatihan, Produktivitas)
## attr(,"factors")
##               Pelatihan Produktivitas Pelatihan:Produktivitas
## Pelatihan             1             0                       1
## Produktivitas         0             1                       1
## attr(,"term.labels")
## [1] "Pelatihan"               "Produktivitas"          
## [3] "Pelatihan:Produktivitas"
## 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
Statistik Makna
X² = 0 Tidak ada deviasi/residual, artinya model pas dengan data.
df = 0 Karena ini model saturated, tidak ada derajat kebebasan tersisa. Model sepenuhnya “menghafal” data.
P-value = 1 P-value maksimum menunjukkan tidak ada sisa kesalahan — model sangat sesuai dengan data.

14.3 Model Independent

Model independen adalah model log-linear yang mengasumsikan tidak ada interaksi antara variabel kategorik. Dengan kata lain, variabel-variabel pada model saling independen satu sama lain.

Misalkan ada dua variabel kategorik \(A\) dan \(B\), model independen dinyatakan sebagai:

\(\log(m_{ij}) = \mu + \lambda_i^A + \lambda_j^B\)

di mana:

\(m_{ij}\) = nilai ekspektasi frekuensi pada sel ke- \((i,j)\)

\(\mu\)= intercept (konstanta)

\(\lambda_i^A\) = efek dari kategori ke- \(i\) pada variabel \(A\)

\(\lambda_j^B\)​ = efek dari kategori ke- \(j\) pada variabel \(B\)

model_ind <- loglm(~ Pelatihan+Produktivitas, data = data_pelatihan)
summary(model_ind)
## Formula:
## ~Pelatihan + Produktivitas
## attr(,"variables")
## list(Pelatihan, Produktivitas)
## attr(,"factors")
##               Pelatihan Produktivitas
## Pelatihan             1             0
## Produktivitas         0             1
## attr(,"term.labels")
## [1] "Pelatihan"     "Produktivitas"
## 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 25.67248  1 4.04555e-07
## Pearson          25.06394  1 5.54604e-07

Nilai p-value yang didapatkan menunjukkan bahwa data tidak cocok dengan model independen (adanya interaksi).

14.4 Odds Ratio, Log Odds, dan Interpretasi

Nilai Odds Ratio dihasilkan dari perhitungan :

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

Atau dengan kata lain, odds ratio didapatkan dari perbandingan antara hasil perkalian diagonal utama ( 𝑎 dan 𝑑 ) dengan hasil perkalian diagonal sekunder ( 𝑏 dan 𝑐 ).

Interpretasi nilai OR:

• OR = 1: Tidak ada asosiasi

• OR > 1: Asosiasi positif

• OR < 1: Asosiasi negatif

Perhitungan (Dari Data Pelatihan)

OR <- (data_pelatihan[1,1] * data_pelatihan[2,2]) / (data_pelatihan[1,2] * data_pelatihan[2,1])
OR
## [1] 4.5

Dari hasil perhitungan, terlihat bahwa odds ratio sebesar 4.5 menunjukkan bahwa Peluang seseorang yang ikut pelatihan untuk memiliki produktivitas tinggi adalah sekitar 4.5 kali lebih besar dibandingkan dengan yang tidak ikut pelatihan.

Log Odds

logOR <- log(OR)
logOR
## [1] 1.504077

Log Odds Ratio (logOR) menunjukkan arah dan kekuatan hubungan antara dua variabel.

Jika:

  • logOR > 0 : hubungan positif (kejadian lebih mungkin terjadi di kelompok pertama)

  • logOR < 0 : hubungan negatif

  • logOR = 0 : tidak ada perbedaan (OR = 1)

Dari output, dapat diinterpretasikan bahwa terdapat hubungan positif yang kuat antara dua variabel (Produktivitas & Pelatihan).

Confidence Interval

# Hitung standard error log(OR)
SE_logOR <- sqrt(1/data_pelatihan[1,1] + 1/data_pelatihan[1,2] + 1/data_pelatihan[2,1] + 1/data_pelatihan[2,2])

# Hitung CI log(OR)
z <- 1.96  # untuk CI 95%
lower_log <- logOR - z * SE_logOR
upper_log <- logOR + z * SE_logOR

# Transformasi balik ke OR
lower_CI <- exp(lower_log)
upper_CI <- exp(upper_log)

# Output hasil
cat("Odds Ratio =", OR, "\n")
## Odds Ratio = 4.5
cat("95% CI for OR = [", lower_CI, ",", upper_CI, "]\n")
## 95% CI for OR = [ 2.459517 , 8.233325 ]

Interpretasi : Terdapat keyakinan sebesar 95% bahwa nilai odds ratio populasi terletak antara batas bawah dan atas dari interval ini.

14.5 Perbandingan Model

anova(model_ind, model_saturated)
## LR tests for hierarchical log-linear models
## 
## Model 1:
##  ~Pelatihan + Produktivitas 
## Model 2:
##  ~Pelatihan * Produktivitas 
## 
##           Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1   25.67248  1                                    
## Model 2    0.00000  0   25.67248         1              0
## Saturated  0.00000  0    0.00000         0              1

Interaksi antara Pelatihan dan Produktivitas signifikan secara statistik (p < 0.001), sehingga model dengan interaksi (model saturated) lebih baik digunakan daripada model tanpa interaksi (model independen).

15 Model Log Linear Tiga Arah

Model log-linear tiga arah merupakan salah satu tipe model statistik yang digunakan untuk menyelidiki keterkaitan antara tiga variabel kategorik dalam tabel kontingensi dengan tiga dimensi. Model ini sangat bermanfaat untuk menemukan hubungan serta interaksi antara ketiga variabel tersebut.

15.1 Model Log-linear untuk Tabel Tiga Arah

15.1.1 Model Saturated

\[ \log (\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} + \lambda_{ijk}^{XYZ} \]

Model ini mencakup semua efek utama yang mungkin terjadi, interaksi dua arah, dan interaksi tiga arah, yang secara sempurna sesuai dengan frekuensi yang diamati.

15.1.2 Model Homogenous

\[ \log (\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

Model ini mencakup semua efek utama dan interaksi dua arah, tetapi meniadakan interaksi tiga arah variabel \(X\), \(Y\), \(Z\) \((\lambda_{ijk}^{XYZ})\).

15.1.3 Model Conditional

Conditional on \(X\)

\[ \log (\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} \]

Model ini mengasumsikan bahwa \(Y\) dan \(Z\) adalah independen bersyarat terhadap \(X\), termasuk semua efek utama dan interaksi yang melibatkan \(X\).

Conditional on \(Y\)

\[ \log (\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{jk}^{YZ} \]

Model ini mengasumsikan bahwa \(X\) dan \(Z\) saling bebas (tidak berhubungan) jika diketahui \(Y\), termasuk efek utama dan interaksi yang melibatkan \(Y\).

Conditional on \(Z\)

\[ \log (\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

Model ini mengasumsikan bahwa \(X\) dan \(Y\) saling bebas (tidak berhubungan) jika diketahui \(Z\), termasuk efek utama dan interaksi yang melibatkan \(Z\)

15.1.4 Model Joint Independence

Independensi \(X\) dan \(Y\)

\[ \log (\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} \]

Model ini mengasumsikan bahwa pasangan variabel \(X\) dan \(Y\) secara bersama-sama tidak berhubungan dengan \(Z\). Model tetap memasukkan efek utama serta interaksi antara \(X\) dan \(Y\), tetapi tidak ada interaksi yang melibatkan \(Z\).

Independensi \(X\) dan \(Z\)

\[ \log (\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ik}^{XZ} \]

Model joint independence antara \(X\) dan \(Z\) terhadap \(Y\) menyatakan bahwa kombinasi dari variabel \(X\) dan \(Z\) secara bersama-sama tidak berhubungan dengan \(Y\).

Independensi \(Y\) dan \(Z\)

\[ \log (\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{jk}^{YZ} \]

Model joint independence antara \(Y\) dan \(Z\) terhadap \(X\) menyatakan bahwa kombinasi dari variabel \(Y\) dan \(Z\) secara bersama-sama tidak berhubungan dengan \(X\).

Model Tanpa Interaksi

\[ \log(\mu_{ijk}) = \mu + \lambda_i^X + \lambda_j^Y + \lambda_k^Z \]

Model log-linear tanpa interaksi (saling bebas) mengasumsikan semua variabel (\(X\), \(Y\), \(Z\)) adalah independen. Model ini hanya mencakup efek utama, dengan asumsi tidak ada interaksi antara variabel

15.2 Pengujian Interaksi

Dalam analisis model log-linear tiga arah, pengujian interaksi bertujuan untuk mendeteksi keberadaan atau ketiadaan interaksi di antara variabel-variabel yang terlibat. Proses pengujian ini dilakukan secara berurutan, dimulai dari interaksi tingkat tertinggi hingga ke tingkat yang lebih rendah. Untuk model log-linear yang melibatkan tiga variabel, misalnya X, Y, dan Z, tahapan pengujian mencakup langkah-langkah berikut:

15.2.1 Pengujian Interaksi Tiga Arah (XYZ):

  • Perbandingan Model: Membandingkan model saturated (yang mencakup semua interaksi) dengan model homogenous (yang mengasumsikan keseragaman tertentu dalam interaksi).

15.2.2 Pengujian Interaksi Dua Arah (XY, XZ, YZ):

  • Langkah 1: Membandingkan model homogenous dengan model conditional, untuk mengevaluasi keberadaan interaksi dua arah tertentu.
  • Langkah 2: Membandingkan model conditional dengan model joint independence, yang menguji independensi bersama antar variabel.
  • Langkah 3: Membandingkan model joint independence dengan model tanpa interaksi, yang mengasumsikan tidak adanya hubungan interaktif antar variabel.

Setiap langkah dalam proses pengujian ini dilakukan untuk menilai sejauh mana model sesuai dengan data yang ada serta untuk mengidentifikasi struktur interaksi yang paling tepat mencerminkan pola dalam data yang diamati. Dengan pendekatan ini, analisis bertujuan menemukan model yang paling cocok untuk menjelaskan hubungan antar variabel berdasarkan observasi.

15.3 Analisis Log-linear Untuk Tabel Tiga Arah

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

Simulasi Data

Data pada contoh kasus ini merupakan simulasi kontingensi 3 arah dengan kategori Jenis Kelamin, Usia, dan Hiburan.

jeniskelamin <- c("Pria", "Wanita")
usia <- c("Pemuda", "Dewasa", "Lansia")
hiburan <- c("Film", "Musik", "Olahraga")
frequencies <- c(
  60, 100,  
  50, 70,   
  30, 50,
  80, 50,   
  60, 70,   
  40, 80,   
  80, 40,   
  60, 60,   
  40, 60    
  )
            
cont_table <- array(frequencies, dim = c(3, 2, 3), 
                    dimnames = list(Usia = usia, JenisKelamin = jeniskelamin, Hiburan = hiburan))
print(cont_table)
## , , Hiburan = Film
## 
##         JenisKelamin
## Usia     Pria Wanita
##   Pemuda   60     70
##   Dewasa  100     30
##   Lansia   50     50
## 
## , , Hiburan = Musik
## 
##         JenisKelamin
## Usia     Pria Wanita
##   Pemuda   80     70
##   Dewasa   50     40
##   Lansia   60     80
## 
## , , Hiburan = Olahraga
## 
##         JenisKelamin
## Usia     Pria Wanita
##   Pemuda   80     60
##   Dewasa   40     40
##   Lansia   60     60
ftable(cont_table, row.vars = c("Usia", "JenisKelamin"), col.vars = "Hiburan")
##                     Hiburan Film Musik Olahraga
## Usia   JenisKelamin                            
## Pemuda Pria                   60    80       80
##        Wanita                 70    70       60
## Dewasa Pria                  100    50       40
##        Wanita                 30    40       40
## Lansia Pria                   50    60       60
##        Wanita                 50    80       60

Penentuan Referensi

Dalam model log-linear, kategori referensi (atau reference category) diperlukan untuk mendefinisikan parameter dalam model agar dapat diestimasi secara unik (mengatasi masalah overparameterization)

Referensi yang digunakan :

  • Jenis Kelamin: Jika Pria dipilih sebagai referensi, efek Wanita akan diestimasi relatif terhadap Pria. Ini masuk akal jika Pria dianggap kelompok mayoritas atau baseline.
  • Usia: Memilih Pemuda sebagai referensi memungkinkan perbandingan Dewasa dan Lansia terhadap kelompok termuda, yang relevan untuk studi perkembangan preferensi.
  • Hiburan: Memilih Film sebagai referensi masuk akal jika Film adalah hiburan paling umum atau standar untuk perbandingan dengan Musik dan Olahraga.
jeniskelamin <- relevel(factor(jeniskelamin), ref = "Pria")
usia <- relevel(factor(usia), ref = "Pemuda")
hiburan <- relevel(factor(hiburan), ref = "Film")

15.3.1 Model Saturated

Model saturated memuat semua interaksi baik dua maupun tiga arah (dalam kasus kontingensi tiga arah).

data_hiburan <- expand.grid(Usia = usia, JK = jeniskelamin, Hib = hiburan)
data_hiburan$Freq = as.vector(cont_table)
model_saturated <- glm(Freq ~ Usia * JK * Hib, 
                       family = poisson(link = "log"), 
                       data = data_hiburan)
summary(model_saturated)
## 
## Call:
## glm(formula = Freq ~ Usia * JK * Hib, family = poisson(link = "log"), 
##     data = data_hiburan)
## 
## Coefficients:
##                                 Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       4.0943     0.1291  31.715  < 2e-16 ***
## UsiaDewasa                        0.5108     0.1633   3.128 0.001759 ** 
## UsiaLansia                       -0.1823     0.1915  -0.952 0.341024    
## JKWanita                          0.1542     0.1759   0.876 0.380927    
## HibMusik                          0.2877     0.1708   1.684 0.092086 .  
## HibOlahraga                       0.2877     0.1708   1.684 0.092086 .  
## UsiaDewasa:JKWanita              -1.3581     0.2726  -4.983 6.26e-07 ***
## UsiaLansia:JKWanita              -0.1542     0.2664  -0.579 0.562784    
## UsiaDewasa:HibMusik              -0.9808     0.2432  -4.032 5.52e-05 ***
## UsiaLansia:HibMusik              -0.1054     0.2566  -0.411 0.681341    
## UsiaDewasa:HibOlahraga           -1.2040     0.2533  -4.753 2.00e-06 ***
## UsiaLansia:HibOlahraga           -0.1054     0.2566  -0.411 0.681341    
## JKWanita:HibMusik                -0.2877     0.2403  -1.197 0.231213    
## JKWanita:HibOlahraga             -0.4418     0.2452  -1.802 0.071547 .  
## UsiaDewasa:JKWanita:HibMusik      1.2685     0.3822   3.319 0.000903 ***
## UsiaLansia:JKWanita:HibMusik      0.5754     0.3562   1.615 0.106285    
## UsiaDewasa:JKWanita:HibOlahraga   1.6458     0.3917   4.201 2.65e-05 ***
## UsiaLansia:JKWanita:HibOlahraga   0.4418     0.3653   1.209 0.226483    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 8.9818e+01  on 17  degrees of freedom
## Residual deviance: 4.0856e-14  on  0  degrees of freedom
## AIC: 142.06
## 
## Number of Fisher Scoring iterations: 3

Mayoritas efek interaksi dua arah maupun tiga arah adalah signifikan (p < 0.05). Oleh karena itu model saturated perlu dipertimbangkan dalam penentuan dan perbandingan model pada tahap selanjutnya.

Koefisien

exp(model_saturated$coefficients)
##                     (Intercept)                      UsiaDewasa 
##                      60.0000000                       1.6666667 
##                      UsiaLansia                        JKWanita 
##                       0.8333333                       1.1666667 
##                        HibMusik                     HibOlahraga 
##                       1.3333333                       1.3333333 
##             UsiaDewasa:JKWanita             UsiaLansia:JKWanita 
##                       0.2571429                       0.8571429 
##             UsiaDewasa:HibMusik             UsiaLansia:HibMusik 
##                       0.3750000                       0.9000000 
##          UsiaDewasa:HibOlahraga          UsiaLansia:HibOlahraga 
##                       0.3000000                       0.9000000 
##               JKWanita:HibMusik            JKWanita:HibOlahraga 
##                       0.7500000                       0.6428571 
##    UsiaDewasa:JKWanita:HibMusik    UsiaLansia:JKWanita:HibMusik 
##                       3.5555556                       1.7777778 
## UsiaDewasa:JKWanita:HibOlahraga UsiaLansia:JKWanita:HibOlahraga 
##                       5.1851852                       1.5555556

Interpretasi koefisien model log-linear saturated yang dieksponensialkan memberikan rasio expected counts (atau multiplicative effect) relatif terhadap kategori referensi. Dalam hal ini, kategori referensi adalah Usia: Pemuda, Jenis Kelamin: Pria, dan Hiburan: Film. Koefisien yang dieksponensialkan mencerminkan bagaimana expected counts di sel tertentu berubah dibandingkan dengan sel referensi.

Intercept: 60.0000000

  • Expected count untuk kombinasi referensi (Pemuda, Pria, Film) adalah 60. Ini adalah frekuensi yang diharapkan di sel referensi ketika semua variabel berada pada kategori baseline.

UsiaDewasa: 1.6666667

  • Untuk Pria dan Film, expected count untuk Dewasa adalah 1.66 kali lipat dari expected count untuk Pemuda.

UsiaDewasa:JKWanita: 0.2571

  • Rasio expected count Wanita/Pria untuk Dewasa, Film adalah 0.257 kali rasio Wanita/Pria untuk Pemuda, Film. Ini menunjukkan penurunan relatif expected count untuk Wanita pada kelompok Dewasa.

UsiaDewasa:HibMusik: 0.375

  • Untuk Pria dan Dewasa, expected count untuk Musik adalah 0.375 kali lipat dari expected count untuk Film dibandingkan dengan Pemuda, Pria, Musik relatif terhadap Pemuda, Pria, Film. Dengan kata lain, expected count untuk Dewasa, Pria, Musik relatif terhadap Dewasa, Pria, Film lebih rendah dibandingkan rasio serupa untuk Pemuda.

UsiaDewasa:HibOlahraga: 0.30

  • Untuk Pria dan Dewasa, expected count untuk Olahraga adalah 0.3 kali lipat dari expected count untuk Film dibandingkan dengan Pemuda, Pria, Olahraga relatif terhadap Pemuda, Pria, Film. Ini menunjukkan penurunan relatif expected count untuk Olahraga pada kelompok Dewasa.

UsiaDewasa:JKWanita:HibMusik: 3.5555556

  • Untuk Dewasa, rasio expected count Wanita/Pria untuk Musik relatif terhadap Film adalah 3.55 kali lipat dari rasio yang sama untuk Pemuda. Dengan kata lain, interaksi tiga arah ini meningkatkan expected count untuk Dewasa, Wanita, Musik secara signifikan dibandingkan dengan kombinasi referensi.

UsiaDewasa:JKWanita:HibOlahraga: 5.1851852

  • Untuk Dewasa, rasio expected count Wanita/Pria untuk Olahraga relatif terhadap Film adalah 5.1851852 kali lipat dari rasio yang sama untuk Pemuda. Ini menunjukkan peningkatan besar dalam expected count untuk Dewasa, Wanita, Olahraga dibandingkan kombinasi referensi, mencerminkan efek interaksi tiga arah yang kuat.

Kesimpulan

Terlihat dari uji signifikansi interaksi bahwa baik model saturated sangat penting dipertimbangkan mengingat banyak dari interaksi yang signifikan (p < 0.05).

15.3.2 Model Homogenous

Dalam kasus kontingensi tiga arah, model homogenus hanya memfokuskan interaksi dua variabel dan mengabaikan interaksi ketiga variabel. Sehingga interaksi tersebut tidak muncul pada model ini

# Jalankan model log-linear homogenous
model_homogenous <- glm(Freq ~ Usia + JK + Hib + 
                        Usia:JK + Usia:Hib + JK:Hib, 
                        family = poisson(link = "log"), 
                        data = data_hiburan)

# Tampilkan ringkasan model
summary(model_homogenous)
## 
## Call:
## glm(formula = Freq ~ Usia + JK + Hib + Usia:JK + Usia:Hib + JK:Hib, 
##     family = poisson(link = "log"), data = data_hiburan)
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)             4.288703   0.105674  40.584  < 2e-16 ***
## UsiaDewasa              0.165725   0.138158   1.200  0.23032    
## UsiaLansia             -0.354506   0.149172  -2.376  0.01748 *  
## JKWanita               -0.243408   0.134138  -1.815  0.06958 .  
## HibMusik                0.015150   0.139301   0.109  0.91340    
## HibOlahraga             0.003723   0.140902   0.026  0.97892    
## UsiaDewasa:JKWanita    -0.426989   0.155556  -2.745  0.00605 ** 
## UsiaLansia:JKWanita     0.198555   0.144117   1.378  0.16828    
## UsiaDewasa:HibMusik    -0.482823   0.182974  -2.639  0.00832 ** 
## UsiaLansia:HibMusik     0.179976   0.177896   1.012  0.31168    
## UsiaDewasa:HibOlahraga -0.543866   0.187869  -2.895  0.00379 ** 
## UsiaLansia:HibOlahraga  0.100611   0.182362   0.552  0.58115    
## JKWanita:HibMusik       0.270494   0.150111   1.802  0.07155 .  
## JKWanita:HibOlahraga    0.153526   0.154408   0.994  0.32008    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 89.818  on 17  degrees of freedom
## Residual deviance: 21.131  on  4  degrees of freedom
## AIC: 155.19
## 
## Number of Fisher Scoring iterations: 4

Dari efek - efek interaksi dua arah pada model ini terlihat bahwa hanya beberapa interaksi yang signifikan. Hal ini bisa disebabkan karena model saturated (model yang memuat semua interaksi dua dan tiga arah) lebih banyak.

15.3.3 Pengujian Hipotesis Interaksi Tiga Arah (Saturated vs Homogenous)

Hipotesis :

H0: Tidak ada interaksi tiga arah (model homogenous lebih cocok)

H1: Ada interaksi tiga arah (model saturated diperlukan)

Statistik Uji :

Deviance dan Chi-square :

anova(model_homogenous, model_saturated, test = "Chisq")

Keputusan :

Tolak H0 (p < 0.05)

Kesimpulan :

Terdapat interaksi tiga arah (antara Usia (Kelompok Umur), jenis kelamin, , dan Preferensi Hiburan) atau dengan kata lain model saturated diperlukan.

Model Conditional on X (Usia)

Model hanya menyertakan interaksi yang melibatkan Usia, yaitu Usia:JenisKelamin (XY) dan Usia:Hiburan (XZ), mengasumsikan Jenis Kelamin dan Hiburan independen kondisional pada Usia (tidak ada interaksi YZ).

model_cond_usia <- glm(Freq ~ Usia + JK + Hib + 
                       Usia:JK + Usia:Hib, 
                       family = poisson(link = "log"), data = data_hiburan)
summary(model_cond_usia)
## 
## Call:
## glm(formula = Freq ~ Usia + JK + Hib + Usia:JK + Usia:Hib, family = poisson(link = "log"), 
##     data = data_hiburan)
## 
## Coefficients:
##                        Estimate Std. Error z value Pr(>|z|)    
## (Intercept)             4.22091    0.09928  42.515  < 2e-16 ***
## UsiaDewasa              0.18987    0.13957   1.360  0.17370    
## UsiaLansia             -0.36604    0.15153  -2.416  0.01571 *  
## JKWanita               -0.09531    0.09770  -0.976  0.32930    
## HibMusik                0.14310    0.11983   1.194  0.23240    
## HibOlahraga             0.07411    0.12180   0.608  0.54289    
## UsiaDewasa:JKWanita    -0.45123    0.15459  -2.919  0.00351 ** 
## UsiaLansia:JKWanita     0.20654    0.14384   1.436  0.15105    
## UsiaDewasa:HibMusik    -0.51083    0.18211  -2.805  0.00503 ** 
## UsiaLansia:HibMusik     0.19337    0.17749   1.089  0.27594    
## UsiaDewasa:HibOlahraga -0.55962    0.18716  -2.990  0.00279 ** 
## UsiaLansia:HibOlahraga  0.10821    0.18212   0.594  0.55239    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 89.818  on 17  degrees of freedom
## Residual deviance: 24.392  on  6  degrees of freedom
## AIC: 154.45
## 
## Number of Fisher Scoring iterations: 4

15.3.4 Pengujian Interaksi antara Y (Usia) dan Z (Hiburan)

(Homogenous vs Conditional on X)

Hipotesis :

H0 : Tidak ada interaksi antara Usia (Y) dan Hiburan (Z)

H1 : Terdapat interaksi antara Usia (Y) dan Hiburan (Z)

Statistik Uji :

Deviance dan Chi-square :

anova(model_cond_usia, model_homogenous, test = "Chisq")

Keputusan :

Terima H0 (p > 0.05)

Kesimpulan :

Tidak ada bukti statistik yang cukup untuk menolak hipotesis nol bahwa interaksi JenisKelamin:Hiburan (JK:Hib) tidak signifikan. Dengan kata lain, menambahkan interaksi JK:Hib ke dalam model (dari conditional on Usia ke homogenous) tidak memberikan perbaikan kecocokan model yang signifikan secara statistik. Ini menunjukkan bahwa hubungan antara Jenis Kelamin dan Hiburan tidak signifikan setelah mengontrol Usia.

Model Conditional on Y (Jenis Kelamin)

Model menyertakan Usia:JenisKelamin (XY) dan JenisKelamin:Hiburan (YZ), mengasumsikan Usia dan Hiburan independen kondisional pada Jenis Kelamin (tidak ada XZ).

model_cond_jk <- glm(Freq ~ Usia + JK + Hib + 
                     Usia:JK + JK:Hib, 
                     family = poisson(link = "log"), data = data_hiburan)
summary(model_cond_jk)
## 
## Call:
## glm(formula = Freq ~ Usia + JK + Hib + Usia:JK + JK:Hib, family = poisson(link = "log"), 
##     data = data_hiburan)
## 
## Coefficients:
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)           4.37771    0.08708  50.271  < 2e-16 ***
## UsiaDewasa           -0.14660    0.09904  -1.480  0.13880    
## UsiaLansia           -0.25783    0.10212  -2.525  0.01157 *  
## JKWanita             -0.28336    0.13134  -2.157  0.03097 *  
## HibMusik             -0.10008    0.10013  -1.000  0.31751    
## HibOlahraga          -0.15415    0.10157  -1.518  0.12911    
## UsiaDewasa:JKWanita  -0.45123    0.15459  -2.919  0.00351 ** 
## UsiaLansia:JKWanita   0.20654    0.14384   1.436  0.15105    
## JKWanita:HibMusik     0.33647    0.14817   2.271  0.02316 *  
## JKWanita:HibOlahraga  0.21869    0.15243   1.435  0.15137    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 89.818  on 17  degrees of freedom
## Residual deviance: 38.550  on  8  degrees of freedom
## AIC: 164.61
## 
## Number of Fisher Scoring iterations: 4

15.3.5 Pengujian Interaksi antara X (Jenis Kelamin) dan Z (Hiburan)

(Homogenous vs Conditional on Y)

Hipotesis :

H0 : Tidak ada interaksi antara X (Jenis Kelamin) dan Z (Hiburan)

H1 : Terdapat interaksi antara X (Jenis Kelamin) dan Z (Hiburan)

Statistik Uji :

Deviance dan Chi-square

anova(model_cond_jk, model_homogenous, test = "Chisq")

Keputusan :

Tolak H0 (p < 0.05)

Kesimpulan :

Penambahan interaksi dua arah Usia:Hiburan dalam model homogenous secara signifikan meningkatkan kecocokan model dibandingkan model conditional on Jenis Kelamin. Menghilangkan interaksi Usia:Hiburan menyebabkan kehilangan informasi signifikan tentang hubungan antar variabel, sehingga model homogenous (yang menyertakan interaksi Usia:Hiburan) lebih sesuai untuk data ini.

Model Conditional on Z (Hiburan)

Model menyertakan Usia:Hiburan (XZ) dan JenisKelamin:Hiburan (YZ), mengasumsikan Usia dan Jenis Kelamin independen kondisional pada Hiburan (tidak ada XY).

#Model conditional on Hiburan (tanpa Usia:JenisKelamin)
model_cond_hib <- glm(Freq ~ Usia + JK + Hib + 
                      Usia:Hib + JK:Hib, 
                      family = poisson(link = "log"), data = data_hiburan)
summary(model_cond_hib)
## 
## Call:
## glm(formula = Freq ~ Usia + JK + Hib + Usia:Hib + JK:Hib, family = poisson(link = "log"), 
##     data = data_hiburan)
## 
## Coefficients:
##                          Estimate Std. Error z value Pr(>|z|)    
## (Intercept)             4.329e+00  9.837e-02  44.003  < 2e-16 ***
## UsiaDewasa              3.450e-11  1.240e-01   0.000  1.00000    
## UsiaLansia             -2.624e-01  1.330e-01  -1.972  0.04855 *  
## JKWanita               -3.365e-01  1.069e-01  -3.147  0.00165 ** 
## HibMusik               -1.105e-02  1.377e-01  -0.080  0.93606    
## HibOlahraga            -2.288e-02  1.394e-01  -0.164  0.86961    
## UsiaDewasa:HibMusik    -5.108e-01  1.821e-01  -2.805  0.00503 ** 
## UsiaLansia:HibMusik     1.934e-01  1.775e-01   1.089  0.27594    
## UsiaDewasa:HibOlahraga -5.596e-01  1.872e-01  -2.990  0.00279 ** 
## UsiaLansia:HibOlahraga  1.082e-01  1.821e-01   0.594  0.55239    
## JKWanita:HibMusik       3.365e-01  1.482e-01   2.271  0.02316 *  
## JKWanita:HibOlahraga    2.187e-01  1.524e-01   1.435  0.15137    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 89.818  on 17  degrees of freedom
## Residual deviance: 36.826  on  6  degrees of freedom
## AIC: 166.89
## 
## Number of Fisher Scoring iterations: 4

15.3.5 Pengujian Interaksi antara X (Jenis Kelamin) dan Y (Usia)

(Homogenous vs Conditional on Z)

Hipotesis :

H0 : Tidak ada interaksi antara X (Jenis Kelamin) dan Y (Usia)

H1 : Terdapat interaksi antara X (Jenis Kelamin) dan Y (Usia)

Statistik Uji :

Deviance dan Chi-square

anova(model_cond_hib, model_homogenous, test = "Chisq")

Keputusan :

Tolak H0 (p < 0.05)

Kesimpulan :

Karena p < 0.05 yang artinya interaksi JenisKelamin:Usia signifikan sehingga model yang terbentuk adalah model yang memuat interaksi JenisKelamin:Usia.

15.4 Pemilihan Model Terbaik

Interaksi Pengujian P - value Keputusan Keterangan
XYZ Saturated vs Homogenous 0.0002983 Tolak H0 Interaksi Signifikan
YZ Homogenous vs Conditional X 0.1959 Terima H0 Tidak ada interaksi
XZ Homogenous vs Conditional Y 0.001602 Tolak H0 Interaksi Signifikan
XY Homogenous vs Conditional Z 0.003908 Tolak H0 Interaksi Signifikan

Model terbaik dibentuk dari efek utama dan interaksi yang signifikan yaitu interaksi tiga arah, interaksi usia dan jenis kelamin (XY) , dan interaksi usia dan hiburan (XZ).

model_terbaik <- glm(Freq ~ Usia + JK + Hib + 
                    Usia:JK + Usia:Hib + 
                    Usia:JK:Hib, 
                    family = poisson(link = "log"), data = data_hiburan)
summary(model_terbaik)
## 
## Call:
## glm(formula = Freq ~ Usia + JK + Hib + Usia:JK + Usia:Hib + Usia:JK:Hib, 
##     family = poisson(link = "log"), data = data_hiburan)
## 
## Coefficients:
##                                   Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                      4.094e+00  1.291e-01  31.715  < 2e-16 ***
## UsiaDewasa                       5.108e-01  1.633e-01   3.128 0.001759 ** 
## UsiaLansia                      -1.823e-01  1.915e-01  -0.952 0.341024    
## JKWanita                         1.542e-01  1.759e-01   0.876 0.380927    
## HibMusik                         2.877e-01  1.708e-01   1.684 0.092086 .  
## HibOlahraga                      2.877e-01  1.708e-01   1.684 0.092086 .  
## UsiaDewasa:JKWanita             -1.358e+00  2.726e-01  -4.983 6.26e-07 ***
## UsiaLansia:JKWanita             -1.542e-01  2.664e-01  -0.579 0.562784    
## UsiaDewasa:HibMusik             -9.808e-01  2.432e-01  -4.032 5.52e-05 ***
## UsiaLansia:HibMusik             -1.054e-01  2.566e-01  -0.411 0.681341    
## UsiaDewasa:HibOlahraga          -1.204e+00  2.533e-01  -4.753 2.00e-06 ***
## UsiaLansia:HibOlahraga          -1.054e-01  2.566e-01  -0.411 0.681341    
## UsiaPemuda:JKWanita:HibMusik    -2.877e-01  2.403e-01  -1.197 0.231213    
## UsiaDewasa:JKWanita:HibMusik     9.808e-01  2.972e-01   3.300 0.000966 ***
## UsiaLansia:JKWanita:HibMusik     2.877e-01  2.630e-01   1.094 0.274013    
## UsiaPemuda:JKWanita:HibOlahraga -4.418e-01  2.452e-01  -1.802 0.071547 .  
## UsiaDewasa:JKWanita:HibOlahraga  1.204e+00  3.055e-01   3.941 8.12e-05 ***
## UsiaLansia:JKWanita:HibOlahraga -3.608e-16  2.708e-01   0.000 1.000000    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance:  8.9818e+01  on 17  degrees of freedom
## Residual deviance: -3.9968e-14  on  0  degrees of freedom
## AIC: 142.06
## 
## Number of Fisher Scoring iterations: 3

Koefisien Model

data.frame(koef = model_terbaik$coefficients,
           exp_koef = exp(model_terbaik$coefficients)
           )

Nilai eksponensial dari koefisien menunjukkan peluang seseorang digolongkan kategori tertentu berdasarkan baseline kategorinya.

Term \(e^\beta\) Interpretasi
UsiaDewasa 1.6667 Peluang sesorang dengan Usia = Dewasa ( asumsi Pria & Film diabaikan) 1.67 kali lipat lebih tinggi dibandingkan peluang seseorang merupakan pemuda.
UsiaLansia 0.8333 Peluang sesorang dengan Usia = Lansia ( asumsi Pria & Film diabaikan) adalah 0.83 dibandingkan peluang seseorang merupakan pemuda. Atau, peluang seseroang adalah pemuda yaitu sebesar 1/0.833 = 1.2 kali lebih tinggi dibanding dengan peluang seseorang adalah lansia.
JKWanita 1.1667 Peluang seseorang adalah wanita 1.67 kali lipat lebih tinggi dibandingkan peluang seseorang adalah pria ( asumsi Pemuda & Film diabaikan)
HibMusik 1.3333 Peluang seseorang memilih musik 1.33 kali lipat lebih tinggi dibandingkan peluang seseorang memilih film (asumsi Pemuda & Pria diabaikan)
HibOlahraga 1.3333 Peluang seseorang memilih olahraga 1.33 kali lipat lebih tinggi dibandingkan peluang seseorang memilih film (asumsi Pemuda & Pria diabaikan)

Sedangkan untuk interaksi - interaksi, nilai exponen dari koefisien merupakan nilai odds ratio dibandingkan dengan bsaeline nya dengan memperhatikan asumsi kategori lain diabaikan.

Terlihat bahwa untuk menginterpretasikan model dengan interaksi menimbulkan kompleksitas.

Nilai Dugaan Model Terbaik

data_hiburan$Fitted <- fitted(model_terbaik)

# Tampilkan data frame
print(data_hiburan)
##      Usia     JK      Hib Freq Fitted
## 1  Pemuda   Pria     Film   60     60
## 2  Dewasa   Pria     Film  100    100
## 3  Lansia   Pria     Film   50     50
## 4  Pemuda Wanita     Film   70     70
## 5  Dewasa Wanita     Film   30     30
## 6  Lansia Wanita     Film   50     50
## 7  Pemuda   Pria    Musik   80     80
## 8  Dewasa   Pria    Musik   50     50
## 9  Lansia   Pria    Musik   60     60
## 10 Pemuda Wanita    Musik   70     70
## 11 Dewasa Wanita    Musik   40     40
## 12 Lansia Wanita    Musik   80     80
## 13 Pemuda   Pria Olahraga   80     80
## 14 Dewasa   Pria Olahraga   40     40
## 15 Lansia   Pria Olahraga   60     60
## 16 Pemuda Wanita Olahraga   60     60
## 17 Dewasa Wanita Olahraga   40     40
## 18 Lansia Wanita Olahraga   60     60

Ternyata nilai dugaan persis sama dengan data asli, hal ini kemungkinan terjadi karena model memasukkan interaksi tiga arah, selain itu data simulasi memperkuat risiko overfitting karena frekuensi yang terstruktur (bulat, pola jelas) memudahkan model untuk mencapai kecocokan sempurna. Dalam data nyata dengan lebih banyak variabilitas, model ini mungkin tidak menggeneralisasi dengan baik.

16 Studi Kasus

Studi Kasus 1 : Regresi Logistik

Penggunaan regresi logistik untuk memprediksi probabilitas seorang pasien menderita diabetes berdasarkan usia, Indeks Massa Tubuh (BMI), dan tekanan darah.

Dataset ini terdiri dari 500 pasien dengan variabel-variabel berikut:

Diabetes: Variabel biner (1 jika pasien menderita diabetes, 0 jika tidak)

Age: Usia pasien dalam tahun BMI: Indeks Massa Tubuh dalam kg/m²

BP: Tekanan darah dalam mmHg

Data ini bersifat hipotetis dan disimulasikan untuk keperluan analisis. Spesifikasi Model

Model

Model probabilitas seseorang menderita diabetes dengan regresi logistik :

\[\log\left(\frac{p}{1-p}\right) = \beta_0 + \beta_1 \cdot \text{Age} + \beta_2 \cdot \text{BMI} + \beta_3 \cdot \text{BP}\] di mana (\(p\)) adalah probabilitas pasien menderita diabetes, (\(\beta_0\)) adalah intersep, dan (\(\beta_1\), \(\beta_2\), \(\beta_3\)) adalah koefisien untuk masing-masing variabel prediktor.

Penaksiran Parameter

Parameter dalam regresi logistik ditaksir menggunakan metode maksimum likelihood. Salah satu pendekatan untuk mencapai estimasi ini adalah metode Newton-Raphson, sebuah algoritma iteratif yang memperbarui estimasi parameter dengan rumus: \[\beta^{(t+1)} = \beta^{(t)} - H^{-1}(\beta^{(t)}) \cdot \nabla L(\beta^{(t)})\] di mana (\(H\)) adalah matriks Hessian dari turunan kedua log-likelihood, dan (\(\nabla\)\(L\)) adalah gradien dari fungsi log-likelihood. Dalam regresi logistik, metode ini mengarah pada estimasi parameter yang memaksimalkan likelihood. Dalam praktiknya, akan digunakan fungsi glm di R dengan family = binomial, yang secara internal menerapkan algoritma iteratively reweighted least squares (IRLS), yang sama dengan Newton-Raphson dalam konteks ini.

set.seed(123)
n <- 500
age <- rnorm(n, mean = 48, sd = 14)
bmi <- rnorm(n, mean = 23, sd = 5)
bp <- rnorm(n, mean = 118, sd = 17)

# Linear predictor
eta <- -5 + 0.05 * age + 0.1 * bmi + 0.02 * bp
p <- 1 / (1 + exp(-eta))
diabetes <- rbinom(n, 1, p)

# Data frame
data <- data.frame(diabetes, age, bmi, bp)

# Fit model
model <- glm(diabetes ~ age + bmi + bp, family = binomial, data = data)
summary(model)
## 
## Call:
## glm(formula = diabetes ~ age + bmi + bp, family = binomial, data = data)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -3.366892   1.271426  -2.648  0.00809 ** 
## age          0.027249   0.010481   2.600  0.00933 ** 
## bmi          0.104228   0.029441   3.540  0.00040 ***
## bp           0.015158   0.008464   1.791  0.07332 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 370.89  on 499  degrees of freedom
## Residual deviance: 347.86  on 496  degrees of freedom
## AIC: 355.86
## 
## Number of Fisher Scoring iterations: 5

Uji Wald

# Ambil hasil uji Wald dari model
summary(model)$coefficients
##                Estimate  Std. Error   z value     Pr(>|z|)
## (Intercept) -3.36689158 1.271426107 -2.648122 0.0080940285
## age          0.02724861 0.010481168  2.599769 0.0093286545
## bmi          0.10422848 0.029440843  3.540268 0.0003997203
## bp           0.01515767 0.008464105  1.790818 0.0733225636
wald_test <- summary(model)$coefficients

# Ubah ke dalam data frame untuk tampilan yang lebih baik
wald_test_df <- as.data.frame(wald_test)

# Tambahkan nama variabel sebagai kolom
wald_test_df$Variable <- rownames(wald_test_df)

# Pilih dan urutkan kolom yang relevan
wald_test_df <- wald_test_df[, c("Variable", "Estimate", "Std. Error", "z value", "Pr(>|z|)")]

# Tampilkan tabel
print(wald_test_df)
##                Variable    Estimate  Std. Error   z value     Pr(>|z|)
## (Intercept) (Intercept) -3.36689158 1.271426107 -2.648122 0.0080940285
## age                 age  0.02724861 0.010481168  2.599769 0.0093286545
## bmi                 bmi  0.10422848 0.029440843  3.540268 0.0003997203
## bp                   bp  0.01515767 0.008464105  1.790818 0.0733225636

Koefisien - koefisien yang mempunyai nilai p < 0.05 dianggap sebagai faktor dengan pengaruh yang signifikan terhadap risiko diabetes.

Goodness of Fit

Untuk menilai kecocokan model, akan digunakan uji rasio likelihood (likelihood ratio test) untuk membandingkan model yang disesuaikan dengan model nol (hanya intersep). Selain itu, Akaike Information Criterion (AIC) akan berperan sebagai ukuran kecocokan model.

# Likelihood ratio test
null_model <- glm(diabetes ~ 1, family = binomial, data = data)
anova(null_model, model, test = "Chisq")

Dari hasil yang menunjukkan p < 0.05 yang artinya model lebih baik dari model null.

AIC

# AIC
AIC(model)
## [1] 355.8588

Studi Kasus 2 : Regresi Multinomial

Terdapat data status mahasiswa setelah 6 bulan lulus dari kampus, adapun kategori yang ingin diuji adalah Jenis Kelamin dan Rumpun ilmu dari individu mahasiswa. Akan diselidiki hubungan dan pengaruh Jenis Kelamin dan Rumpun ilmu dengan status mahasiswa setelah 6 bulan lulus.

library(tibble)
## Warning: package 'tibble' was built under R version 4.3.3
status_data <- tribble(
  ~Rumpun,   ~Gender, ~Status,      ~Count,
  "Soshum",  "Wanita",      "Menganggur",      45,
  "Soshum",  "Wanita",      "Kerja",           110,
  "Soshum",  "Wanita",      "Studi Lanjut",    70,
  "Soshum",  "Pria",        "Menganggur",      60,
  "Soshum",  "Pria",        "Kerja",           70,
  "Soshum",  "Pria",        "Studi Lanjut",    30,
  "Saintek", "Wanita",      "Menganggur",      50,
  "Saintek", "Wanita",      "Kerja",           100,
  "Saintek", "Wanita",      "Studi Lanjut",    90,
  "Saintek", "Pria",        "Menganggur",      45,
  "Saintek", "Pria",        "Kerja",           60,
  "Saintek", "Pria",        "Studi Lanjut",    55
)

Ubah menjadi faktor :

status_data <- status_data %>%
  mutate(
    Status = factor(Status, levels = c("Menganggur", "Kerja", "Studi Lanjut")),
    Gender = factor(Gender),
    Rumpun = factor(Rumpun)
  )

Ekspansi Data :

library(tidyr)
expanded_data <- status_data %>%
  uncount(weights = Count)

# Lihat hasil
head(expanded_data)

Model :

library(nnet)
model_mlr <- multinom(Status ~ Rumpun + Gender, data = expanded_data, trace = FALSE)
summary(model_mlr)
## Call:
## multinom(formula = Status ~ Rumpun + Gender, data = expanded_data, 
##     trace = FALSE)
## 
## Coefficients:
##              (Intercept) RumpunSoshum GenderWanita
## Kerja        0.201923614   0.02192909    0.5798228
## Studi Lanjut 0.008092697  -0.46641618    0.7292053
## 
## Std. Errors:
##              (Intercept) RumpunSoshum GenderWanita
## Kerja          0.1622295    0.1797262    0.1802972
## Studi Lanjut   0.1725574    0.1939786    0.1958357
## 
## Residual Deviance: 1661.055 
## AIC: 1673.055

P - value

library(knitr)
library(kableExtra)
## Warning: package 'kableExtra' was built under R version 4.3.3
## 
## Attaching package: 'kableExtra'
## The following object is masked from 'package:dplyr':
## 
##     group_rows
z_values <- summary(model_mlr)$coefficients / summary(model_mlr)$standard.errors
p_values <- 2 * (1 - pnorm(abs(z_values)))

coef_table <- cbind(summary(model_mlr)$coefficients,
                    "z value" = round(z_values, 2),
                    "p value" = round(p_values, 4))
coef_table
##              (Intercept) RumpunSoshum GenderWanita (Intercept) RumpunSoshum
## Kerja        0.201923614   0.02192909    0.5798228        1.24         0.12
## Studi Lanjut 0.008092697  -0.46641618    0.7292053        0.05        -2.40
##              GenderWanita (Intercept) RumpunSoshum GenderWanita
## Kerja                3.22      0.2132       0.9029       0.0013
## Studi Lanjut         3.72      0.9626       0.0162       0.0002
knitr::kable(coef_table, caption = "Koefisien Model Multinomial Logistic Regression")
Koefisien Model Multinomial Logistic Regression
(Intercept) RumpunSoshum GenderWanita (Intercept) RumpunSoshum GenderWanita (Intercept) RumpunSoshum GenderWanita
Kerja 0.2019236 0.0219291 0.5798228 1.24 0.12 3.22 0.2132 0.9029 0.0013
Studi Lanjut 0.0080927 -0.4664162 0.7292053 0.05 -2.40 3.72 0.9626 0.0162 0.0002

Prediksi Peluang

expanded_data$Rumpun <- factor(expanded_data$Rumpun)
expanded_data$Gender <- factor(expanded_data$Gender)

new_data <- expand.grid(
  Rumpun = levels(expanded_data$Rumpun),
  Gender = levels(expanded_data$Gender)
)
# Prediksi probabilitas
pred_probs <- predict(model_mlr, newdata = new_data, type = "probs")

# Gabungkan
final_result <- cbind(new_data, pred_probs)
print(final_result)
##    Rumpun Gender Menganggur     Kerja Studi Lanjut
## 1 Saintek   Pria  0.3094174 0.3786510    0.3119316
## 2  Soshum   Pria  0.3468333 0.4338492    0.2193176
## 3 Saintek Wanita  0.1895532 0.4142277    0.3962191
## 4  Soshum Wanita  0.2200290 0.4914866    0.2884844

Dapat dilihat bahwa status dengan peluang terbesar pada masing - masing kategori adalah status bekerja.

# Model null: hanya intercept
model_null <- multinom(Status ~ 1, data = expanded_data, trace = FALSE)

# Log-likelihood
ll_null <- logLik(model_null)
ll_full <- logLik(model_mlr)

# G2 (deviance)
G2 <- -2 * (as.numeric(ll_null) - as.numeric(ll_full))

# Derajat kebebasan
df <- attr(ll_full, "df") - attr(ll_null, "df")

# p-value
pval <- 1 - pchisq(G2, df)

# McFadden's pseudo R2
pseudo_r2 <- 1 - (as.numeric(ll_full) / as.numeric(ll_null))

# Tampilkan hasil
cat("Nilai G² (Deviance):", round(G2, 4), "\n")
## Nilai G² (Deviance): 25.4415

Derajat kebebasan

cat("Derajat Kebebasan:", df, "\n")
## Derajat Kebebasan: 4

p - value

cat("p-value:", round(pval, 4), "\n")
## p-value: 0

Pseudo R

cat("McFadden's Pseudo R²:", round(pseudo_r2, 4), "\n")
## McFadden's Pseudo R²: 0.0151

17 Referensi

Jaya, I.G.N.M, (2025). E-book Analisis Data Kategori.

Hayati, Cut, and Syarifah Nurliana (2022) “Analisis Faktor-Faktor yang Mempengaruhi Minat Mahasiswa dalam Berwirausaha (Studi Kasus pada Mahasiswa Fakultas Ekonomi Universitas Almuslim).” Visioner: Jurnal Ilmiah Ilmu Manajemen dan Kewirausahaan, vol. 14, no. 3, 2022, pp. 171–178. Journal Universitas Almuslim,

Agresti, A. (2013) Categorical Data Analysis. 3rd Edition, John Wiley & Sons Inc., Hoboken.

McCullagh, P., & Nelder, J. A. (1989). Generalized linear models (2nd ed.). Chapman & Hall.​

Fleiss, J.L., Levin, B. and Paik, M.C. (2003) Statistical Methods for Rates and Proportions. 3rd Edition, Wiley, Hoboken. https://doi.org/10.1002/0471445428

Kleinbaum, D. G., & Klein, M. (2010). Logistic Regression. In Statistics for Biology and Health. Springer New York. https://doi.org/10.1007/978-1-4419-1742-3

Christensen, R (1997) Log-Linear Models and Logistic Regression. In Springer Texts in Statistics. Springer-Verlag. https://doi.org/10.1007/b97647

Agresti, A. (2007). An introduction to categorical data analysis. Wiley-Interscience.

Lutur, Jane K., et al (2016). “Perbedaan Pemberian Asi Ekslusif Dan Asi Non Ekslusif Dengan Perubahan Berat Badan Pada Bayi Di Puskesmas Bahu Manado.” Jurnal Keperawatan UNSRAT, vol. 4, no. 2.

Ardiansyahroni, A., Tjalla, A., & Mahdiyah, M. (2023). Data kategorik dalam penelitian: Review bibliometrik. Jurnal Ilmiah Mandala Education9(1).

Astuti, A. B., Efendi, A., Astutik, S., & Sumarminingsih, E. (2020). Analisis Data Kategorik Menggunakan R: Teori dan Aplikasinya pada Berbagai Bidang. Universitas Brawijaya Press.

Sartika, E. (2010). Pengolahan data berskala ordinal. Sigma-Mu2(1), 60-69.

Tamonob, A. M. (2024). ANALISIS DATA KATEGORIK MENGGUNAKAN REGRESI LOGISTIK BINER. REVITALISASI PENGGUNAAN MEDIA SERTA METODE BELAJAR DALAM PEMBELAJARAN MATEMATIKA DAN TEKNIK, 19.