Pertama-tama saya panjatkan puji syukur ke hadirat Tuhan Yang Maha Esa karena atas berkat dan rahmat-Nya, e-book ini dapat saya susun dan selesaikan.
Dalam era serba digital ini, Artificial Intelligent dan Machune Learning bukanlah hal yang baru lagi. Tidak dapat dipungkiri bahwa saat ini kita hidup di era dimana kecerdasan buatan tidak lagi susah untuk ditemui. Salah satu bentuk pemanfaatan AI juga tercermin dalam proses penyusunan e-book ini. Penggunaan AI bukan tanpa sebab, namun sebagai bentuk bahwa AI dapat digunakan sebagai alat bantu dalam kegiatan belajar mengajar.
E-book ini membahas tentang analisis data di bidang data kategori. Terdapat beberapa materi tentang data kategori seperti: Tabel Kontingensi, Regresi Logistik, Ukuran Asosiasi, Model Regresi Logistik, dan Model Log-Linear.
Melalui e-book ini, saya harap banyak mahasiswa maupun praktisi yang terbantu dalam menganalisis data kategori. Akhir kata, saya menyadari bahwa e-book ini memiliki keterbatasan dalam penulisannya, sehingga kritik dan saran akan sangat membantu saya dalam mengembangkan hal yang lebih baik.
Dalam statistik, data kategori adalah jenis data yang dibagi ke dalam kelompok atau kategori dan bukan berupa angka yang bisa diukur. Data seperti ini sering digunakan di berbagai bidang untuk memahami pola, hubungan, dan tren yang tidak bisa dilihat langsung melalui angka. Contoh dari data kategori yaitu jenis kelamin (laki-laki/perempuan), status pernikahan (menikah/belum menikah), tingkat pendidikan (rendah/menengah/tinggi), serta preferensi konsumen (setuju/netral/tidak setuju). Analisis terhadap data kategori sangat penting untuk mengeksplorasi dan memahami informasi yang bersifat nominal atau ordinal.
Karena berbeda dari data numerik, data kategori perlu dianalisis dengan metode khusus, seperti tabel kontingensi, uji chi-square, regresi logistik, dan pendekatan berbasis probabilitas. Seiring berkembangnya teknologi, analisis data kategori juga makin sering digunakan dalam machine learning dan kecerdasan buatan.
Analisis Data Kategori memiliki banyak tujuan, diantaranya yaitu:
Analisis data kategori memungkinkan peneliti untuk melihat bagaimana distribusi data tiap kategori. Misalkan terdapat responden dengan kategori usia muda, dewasa, dan tua. Peneliti dapat melihat sebaran usia dari responden yang ada.
Dalam analisis data kategori, identifikasi pola bertujuan untuk menemukan kecenderungan perbedaan antar kelompok dalam data. Misalnya, melalui visualisasi seperti diagram batang dan diagram lingkaran, kita dapat melihat kecenderungan responden memilih suatu merk dibandingkan dengan merk lain.
Dari tujuan yang telah dijelaskan sebelumnya, peneliti juga dimudahkan dalam pengambilan keputusan karena tren atau polanya sudah diketahui.
Analisis data kategorik merupakan suatu teknik statistik yang digunakan untuk menganalisis variabel yang bersifat kategorik, dimana data yang dianalisis bersifat nominal atau ordinal.
Data yang bersifat nominal tidak memiliki urutan, contohnya seperti jenis kelamin dan status menikah. Sedangkan data yang bersifat ordinal memiliki urutan, contohnya seperti tingkat pendidikan dan preferensi konsumen.
Dalam analisis kategorik, data juga dapat dibedakan berdasarkan jumlah kategorinya. Data biner berarti data hanya memiliki dua kategori, misalnya ya/tidak. Sedangkan data multikategori berarti data tersebut memiliki lebih dari dua kategori.
Data Kategorik | Data Kuantitatif | |
---|---|---|
Bentuk | Kategori | Numerik |
Sifat | Kualitatif | Kuantitatif |
Skala Pengukuran | Skala nominal/ordinal | Skala Interval/ratio |
Analisis data kategori memiliki berbagai manfaat dalam berbagai bidang, diantaranya yaitu:
Bidang | Tujuan | Penerapan |
---|---|---|
Kesehatan | Mengetahui hubungan antara kondisi pasien dan kategori risiko | Melihat hubungan antara pola hidup (Sehat/Tidak Sehat) dengan penderita Diabetes |
Sosial | Menganalisis opini masyarakat | Menganalisis kepuasan terhadap layanan publik (Sangat Puas – Tidak Puas) |
Ekonomi | Mengkategorikan jenis pekerjaan dan status ekonomi | Melihat hubungan antara jenis pekerjaan dan tingkat pengeluaran |
Dalam analisis data kategori, metode yang digunakan bergantung kepada tujuan penelitian. Ada beberapa metode dalam analisis data kategori, diantaranya yaitu:
Uji chi - square bertujuan untuk menguji hubungan antara dua variabel, contohnya yaitu menguji apakah ada hubungan antara jenis kelamin dan preferensi produk.
Uji Fisher’s Exact bertujuan sama seperti uji chi - square, hanya saya uji ini digunakan ketika frekuensi data berjumlah kecil (<5).
Regresi logistik pada analisis data kategorik bertujuan untuk memprediksi kemungkinan suatu kategori terjadi berdasarkan variabel prediktor. Misalnya memprediksi apakah seseorang siap menikah berdasarkan tingkat pendapatan, kesiapan mental dan fisik, serta kematangan emosional.
Tabulasi silang atau cross tabulation memiliki tujuan untuk menampilkan distribusi frekuensi dari dua variabel dalam bentuk tabel. Misalnya menampilkan jumlah laki - laki dan perempuan dalam tiap kelompok penyakit.
Analisis diskriminan bertujuan untuk mengklasifikasikan objek ke dalam kategori berdasarkan beberapa variabel, misalnya menentukan preferensi produk pelanggan jenis kelamin dan pendapatan.
Analisis data kategorik memiliki banyak metode untuk digunakan sesuai tujuan masing-masing. Oleh karena itu, peneliti harus memiliki tujuan dan pemahaman untuk bisa meneliti data kategorik menggunakan metode yang tepat.
Variabel acak kategori adalah variabel yang hanya memiliki beberapa kategori diskrit sebagai hasilnya. Distribusi probabilitas dari variabel ini menggambarkan kemungkinan terjadinya setiap kategori.
Distribusi Bernoulli adalah distribusi probabilitas diskrit yang hanya memiliki dua kemungkinan hasil: sukses (1) atau gagal (0). Distribusi ini digunakan untuk percobaan yang hanya memiliki dua hasil, seperti lemparan koin (kepala atau ekor).
Distribusi Bernoulli memiliki satu parameter, yaitu:
Jika \(X\) adalah peubah acak yang mengikuti distribusi Bernoulli, maka:
\[ X \sim \text{Bernoulli}(p) \]
Probability Mass Function (PMF) dari distribusi Bernoulli didefinisikan sebagai:
\[ P(X = x) = \begin{cases} p & \text{jika } x = 1 \\\\ 1 - p & \text{jika } x = 0 \end{cases} \]
atau bisa ditulis dengan persamaan:
\[ P(X = x) = p^x (1 - p)^{1 - x}, \quad x \in \{0, 1\} \]
set.seed(123) # untuk replikasi hasil
n <- 100
hasil <- rbinom(n, size = 1, prob = 0.8)
# Lihat 10 hasil pertama
head(hasil, 10)
## [1] 1 1 1 0 0 1 1 0 1 1
Distribusi Binomial adalah distribusi probabilitas diskrit yang menggambarkan jumlah keberhasilan dalam sejumlah percobaan Bernoulli yang independen dan identik. Setiap percobaan memiliki dua hasil: sukses atau gagal, dan probabilitas sukses tetap di setiap percobaan.
Jika \(X\) adalah peubah acak yang mengikuti distribusi binomial, maka:
\[ X \sim \text{Binomial}(n, p) \]
dengan:
Probability Mass Functiion (PMF) dari distribusi binomial yaitu:
\[ P(X = k) = \binom{n}{k} p^k (1 - p)^{n - k}, \quad k = 0, 1, 2, ..., n \]
Dengan:
set.seed(123) # agar hasil bisa direproduksi
simulasi <- rbinom(1000, size = 10, prob = 0.9)
head(simulasi, 10)
## [1] 10 8 9 8 7 10 9 8 9 9
Distribusi Multinomial adalah perluasan dari distribusi Binomial yang mamiliki lebih dari dua hasil. Distribusi ini digunakan untuk menghitung probabilitas hasil dari sejumlah percobaan independen dimana tiap percobaan memiliki lebih dari dua kemungkinan hasil.
Jika terdapat \(k\) kategori dan setiap percobaan memiliki probabilitas:
\[ p_1, p_2, ..., p_k \quad \text{dengan} \quad \sum_{i=1}^k p_i = 1 \]
dan dilakukan \(n\) percobaan, maka jumlah hasil dari setiap kategori mengikuti distribusi multinomial sebagai berikut:
\[ (X_1, X_2, ..., X_k) \sim \text{Multinomial}(n; p_1, p_2, ..., p_k) \]
Probability Mass Functiion dari distribusi Multinomial yaitu:
\[ P(X_1 = x_1, ..., X_k = x_k) = \frac{n!}{x_1!x_2!...x_k!} \cdot p_1^{x_1} \cdot p_2^{x_2} \cdots p_k^{x_k} \]
Dengan syarat:
set.seed(123)
multinomial <- rmultinom(n = 1, size = 15, prob = c(0.2, 0.8, 0.6))
multinomial
## [,1]
## [1,] 1
## [2,] 7
## [3,] 7
Distribusi Poisson adalah distribusi probabilitas diskrit yang menggambarkan jumlah kejadian dalam suatu interval waktu tertentu, dengan asumsi bahwa:
Jika \(X\) adalah peubah acak yang mengikuti distribusi Poisson dengan parameter \(\lambda\), maka:
\[ X \sim \text{Poisson}(\lambda) \]
dimana:
Probability Mass Functiion (PMF) dari distribusi Poisson adalah:
\[ P(X = k) = \frac{e^{-\lambda} \lambda^k}{k!}, \quad k = 0, 1, 2, ... \]
Sifat-sifat penting:
set.seed(123)
poisson <- rpois(15, lambda = 5)
poisson
## [1] 4 7 4 8 9 2 5 8 5 5 9 5 6 5 2
Desain sampling adalah tahap awal dalam proses penelitian, terkhusus dalam pengumpulan data. Tujuannya adalah agar data yang dikumpulkan dapat digunakan untuk menarik kesimpulan yang akurat tentang populasi.
Dalam analisis data kategori, desain sampling merujuk pada metode yang digunakan untuk memilih sampel dari populasi agar data yang diperoleh dapat mewakili keseluruhan populasi dengan baik. Desain sampling sangatlah penting karena pemilihan sampel akan memengaruhi validitas hasil penelitian.
Desain sampling dalam analisis data kategori dapat diklasifikasikan menjadi dua pendekatan utama, yaitu prospective sampling dan retrosprective sampling.
Prospective sampling adalah metode pengambilan sampel dimana peneliti mengikuti partisipan dari waktu sekarang ke masa depan untuk mengamati kejadian atau hasil tertentu. Prospective sampling biasanya digunakan dalam studi kohort prospektif. Beberapa jenis desain sampling dalam metode ini meliputi:
Dalam studi eksperimen, peneliti menentukan intervensi, yaitu peserta secara acak dibagi ke dalam kelompok perlakuan (menerima intervensi) dan kelompok kontrol (tidak menerima intervensi atau menerima plasebo). Teknik sampling yang digunakan meliputi:
Studi kohort adalah jenis penelitian dimana peneliti mengamati peserta tanpa intervensi langsung dimana peserta dikelompokkan berdasarkan paparan alami terhadap faktor risiko (misalnya, merokok atau tidak merokok. Jenis sampling yang umum dalam studi yaitu:
Retrospective sampling adalah metode pengambilan sampel dalam penelitian dimana peneliti melihat ke belakang untuk mengidentifikasi subjek berdasarkan outcome yang sudah terjadi, lalu menelusuri paparan atau faktor penyebabnya sebelumnya.
Case-Control Study adalah salah satu jenis studi observasional yang digunakan untuk menyelidiki hubungan antara faktor risiko dan suatu outcome (misalnya, penyakit). Berbeda dengan studi kohort, studi ini dimulai dari outcome dan bekerja mundur untuk menentukan paparan sebelumnya. Teknik sampling yang sering digunakan meliputi:
Dalam studi kohort retrospektif, data historis digunakan untuk mengelompokkan individu berdasarkan paparan kemudian peneliti menganalisis hasil yang terjadi. Teknik sampling yang sering digunakan meliputi :
Tabel kontingensi 2x2 adalah tabel yang digunakan untuk menyajikan data kategorik dari dua variabel yang masing-masing memiliki dua kategori. Tabel ini berguna untuk melihat hubungan atau asosiasi antara dua variabel, misalnya perhitungan risiko relatif (relative risk), perbedaan risiko (risk difference) dan odds ratio.
Struktur umum Tabel Kontingensi 2x2:
Kategori B1 | Kategori B2 | Total | |
---|---|---|---|
Kategori A1 | a | b | a + b |
Kategori A2 | c | d | c + d |
Total | a + c | b + d | n |
Peluang bersama adalah probabilitas kedua variabel terjadi secara bersamaan dalam suatu sel tabel kontingensi:
\[ P(A_i, B_j) = \frac{n_{ij}}{n} \]
Peluang marginal adalah probabilitas suatu variabel tanpa mempertimbangkan variabel lainnya.
\[ P(A_i) = \frac{n_{i.}}{n} \]
\[ P(B_j) = \frac{n_{.j}}{n} \]
Peluang bersyarat adalah probabilitas suatu kejadian terjadi dengan syarat kejadian lain telah terjadi:
\[ P(B_j | A_i) = \frac{P(A_i, B_j)}{P(A_i)} = \frac{n_{ij}}{n_{i.}} \]
Contoh perhitungan: Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.
Penggunaan Sabuk Pengaman | Cedera Fatal | Cedera Tidak Fatal | Total |
---|---|---|---|
Tidak memakai sabuk pengaman | 1,601 | 162,527 | 164,128 |
Memakai sabuk pengaman | 510 | 412,368 | 412,878 |
Total | 2,111 | 574,895 | 577,006 |
Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988
Langkah 1 : Hitung Peluang Bersama
Langkah 2 : Hitung Peluang Marginal
Langkah 3 : Hitung Peluang Bersyarat
\[ P(\text{Fatal} \mid \text{No Seatbelt}) = \frac{1601}{164128} = 0.0098 \]
\[ P(\text{Fatal} \mid \text{Seatbelt}) = \frac{510}{412878} = 0.0012 \]
\[ P(\text{Seatbelt} \mid \text{Fatal}) = \frac{510}{2111} = 0.2416 \]
\[ P(\text{Seatbelt} \mid \text{NonFatal}) = \frac{412368}{574895} = 0.7173 \]
#Data
data <- matrix(c(1601,162527,510,412368), nrow = 2, byrow = TRUE)
colnames(data) <- c("Fatal", "NonFatal")
rownames(data) <- c("No Seatbelt","Seatbelt")
n <- sum(data)
#Peluang Bersama
p_bersama <- data/n
#Peluang Marginal
p_marginal_baris <- rowSums(data)/n
p_marginal_kolumn <- colSums(data)/n
#Peluang Bersyarat
# Peluang bersyarat: P(Fatal | Seatbelt) dan P(Fatal | No Seatbelt)
p_fatal_given_no_seatbelt <- data[1,1] / sum(data[1,])
p_fatal_given_seatbelt <- data[2,1] / sum(data[2,])
# Peluang bersyarat: P(Seatbelt | Fatal) dan P(Seatbelt | NonFatal)
total_fatal <- sum(data[,1])
total_nonfatal <- sum(data[,2])
p_seatbelt_given_fatal <- data[2,1] / total_fatal
p_seatbelt_given_nonfatal <- data[2,2] / total_nonfatal
# Output
p_fatal_given_no_seatbelt
## [1] 0.009754582
p_fatal_given_seatbelt
## [1] 0.001235232
p_seatbelt_given_fatal
## [1] 0.2415917
p_seatbelt_given_nonfatal
## [1] 0.7172927
#Hasil
list(Peluang_Bersama = p_bersama, Peluang_Marginal_Baris = p_marginal_baris, Peluang_Marginal_Kolumn = p_marginal_kolumn)
## $Peluang_Bersama
## Fatal NonFatal
## No Seatbelt 0.002774668 0.2816730
## Seatbelt 0.000883873 0.7146685
##
## $Peluang_Marginal_Baris
## No Seatbelt Seatbelt
## 0.2844476 0.7155524
##
## $Peluang_Marginal_Kolumn
## Fatal NonFatal
## 0.003658541 0.996341459
Intepretasi :
Dikarenakan \[ \mathit{P}(\text{Fatal} \mid \text{No Seatbelt}) > \mathit{P}(\text{Fatal} \mid \text{Seatbelt}) \] maka dapat disimpulkan bahwa penggunaan sabuk pengaman dapat menurunkan risiko fatalitas dalam kecelakaan.
Ukuran asosiasi dalam analisis data kategorik digunakan untuk mengukur kekuatan hubungan antara dua variabel yang bersifat kategorik. Ukuran asosiasi ini penting karena digunakan untuk mengetahui apakah hubungan yang ada signifikan dan seberapa kuat hubungan tersebut. Berikut adalah beberapa ukuran asosiasi:
Risk Difference mengukur selisih risiko antara dua kelompok. Nilai RD menunjukkan perbedaan proporsi kejadian antara kelompok yang terpapar dan tidak terpapar.
Rumus:
\[ RD = \frac{a}{a + b} - \frac{c}{c + d} \] Dengan
Jika 𝐷>0, maka risiko kejadian lebih tinggi di Kelompok 1 dibandingkan Kelompok 2
Jika 𝐷<0, maka risiko kejadian lebih rendah di Kelompok 1 dibandingkan Kelompok 2
Jika \(\mathit{RD} = 0\), maka tidak ada perbedaan risiko antara dua kelompok.
Relative Risk membandingkan kemungkinan terjadinya suatu kejadian pada kelompok yang terpapar terhadap kelompok yang tidak terpapar.
Rumus:
\[ RR = \frac{a / (a + b)}{c / (c + d)} \] Dengan
Jika >1, maka kejadian lebih sering terjadi di Kelompok 1 dibandingkan Kelompok 2.
Jika <1, maka kejadian lebih jarang terjadi di Kelompok 1 dibandingkan Kelompok 2.
Jika \(\mathit{RR} = 1\), maka tidak ada perbedaan risiko antara dua kelompok.
Odds Ratio digunakan untuk membandingkan peluang terjadinya suatu kejadian antara dua kelompok.
Rumus:
\[ OR = \frac{a \cdot d}{b \cdot c} \]
Dengan
Jika 𝑂 >1, maka peluang kejadian lebih besar di Kelompok 1 dibandingkan Kelompok 2.
Jika 𝑂 <1, maka peluang kejadian lebih kecil di Kelompok 1 dibandingkan Kelompok 2.
Jika \(\mathit{OR} = 1\), maka tidak ada perbedaan peluang kejadian antara dua kelompok.
Contoh perhitungan: Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.
Penggunaan Sabuk Pengaman | Cedera Fatal | Cedera Tidak Fatal | Total |
---|---|---|---|
Tidak memakai sabuk pengaman | 1,601 | 162,527 | 164,128 |
Memakai sabuk pengaman | 510 | 412,368 | 412,878 |
Total | 2,111 | 574,895 | 577,006 |
Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988
Risk Difference \[ RD = \frac{a}{a + b} - \frac{c}{c + d} = \frac{1601}{164128} - \frac{510}{412878} = 0.0098 - 0.0012 = 0.0086 \]
Relative Risk \[ \frac{P(\text{Fatal} \mid \text{No Seatbelt})}{P(\text{Fatal} \mid \text{Seatbelt})} = \frac{\frac{1601}{164128}}{\frac{510}{412878}} = \frac{0.009759}{0.001235} \approx 7.897 \]
Odds Ratio \[ OR = \frac{a \cdot d}{b \cdot c} = \frac{1601 \cdot 412368}{162527 \cdot 510} = \frac{660201168}{82988770} = 7.96 \]
# Data
a <- 1601
b <- 162527
c <- 510
d <- 412368
# Risk Difference
RD <- (a / (a + b)) - (c / (c + d))
# Relative Risk
RR <- (a / (a + b)) / (c / (c + d))
# Odds Ratio
OR <- (a * d) / (b * c)
#Hasil
list(Risk_Difference = RD, Relative_Risk = RR, Odds_Ratio = OR)
## $Risk_Difference
## [1] 0.00851935
##
## $Relative_Risk
## [1] 7.896965
##
## $Odds_Ratio
## [1] 7.964905
Interpretasi :
RD = 0.0086 → Perbedaan risiko kecelakaan fatal sebesar 0.86% lebih tinggi pada pengendara yang tidak pakai sabuk.
RR = 7.90 → Pengendara yang tidak memakai sabuk memiliki 8.17 kali risiko lebih besar dalam mengalami kecelakaan fatal.
OR = 7.96 → Peluang kecelakaan fatal 7.96 kali lebih besar pada pengendara yang tidak memakai sabuk.
Inferensi pada statistik merujuk pada proses penarikan kesimpulan tentang suatu populasi berdasarkan data dari sampel. Dalam konteks tabel kontingensi dua arah, inferensi digunakan untuk mengevaluasi dan memahami hubungan antara dua variabel kategorikal yang ditampilkan dalam bentuk tabel matriks.
Tabel kontingensi sendiri menyajikan distribusi frekuensi dari kombinasi kategori dua variabel. Tujuan utama dari inferensi pada tabel kontingensi adalah untuk menganalisis keterkaitan antarvariabel dengan melakukan proses estimasi dan pengujian hipotesis. Secara umum, proses inferensi dalam tabel kontingensi dua arah mencakup dua pendekatan utama yaitu Estimasi dan Pengujian
Estimasi berfungsi untuk memperkirakan parameter populasi berdasarkan sampel. Estimasi dibagi menjadi dua, yaitu:
Estimasi titik digunakan untuk menentukan satu nilai spesifik sebagai perkiraan terbaik dari parameter populasi \[ \hat{p} = \frac{x}{n} \] Dengan:
Estimasi interval berfungsi untuk memberikan rentang nilai yang diyakini mengandung parameter populasi dengan tingkat kepercayaan tertentu
\[ \hat{p} \pm Z_{\alpha/2} \sqrt{ \hat{p}(1 - \hat{p}) } \] Dengan:
Uji ini digunakan untuk membandingkan proporsi kejadian antara dua kelompok dalam tabel kontingensi, terutama untuk menentukan apakah terdapat perbedaan yang signifikan dalam proporsi kejadian antara dua kelompok yang berbeda. Tabel kontingensi 2×2 memiliki struktur sebagai berikut:
Kejadian (+) | Tidak Kejadian (-) | Total (n) | |
---|---|---|---|
Group 1 | n11 | n12 | n1. |
Group 2 | n21 | n22 | n2. |
Total | n.1 | n.2 | n |
Untuk menguji hipotesis bahwa tidak ada perbedaan proporsi antara dua kelompok, kita menggunakan uji z dua proporsi, dengan hipotesis:
- Hipotesis Nol (H0) : tidak ada perbedaan proporsi antara dua kelompok
- Hipotesis Alternatif (H1) : terdapat perbedaan proporsi antara dua kelompok
Estimasi proporsi dalam masing - masing kelompok diberikan oleh :
\[ \hat{p}_1 = \frac{a}{a + b} \]
\[ \hat{p}_2 = \frac{c}{c + d} \]
Estimasi proporsi gabungan
\[ \hat{p} = \frac{x_1 + x_2}{n_1 + n_2} \]
Statistik uji untuk uji proporsi dua sampel :
\[ Z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{\hat{p}(1 - \hat{p})\left(\frac{1}{n_1} + \frac{1}{n_2}\right)}} \]
Setelah nilai Z diketahui, kita bisa melihat keputusan apakah tolak H0 atau terima H0 dengan cara apabila |Z| lebih besar dari nilai kritis tertentu untuk tingkat signifikansi tertentu, maka H0 ditolak.
Contoh perhitungan: Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.
Penggunaan Sabuk Pengaman | Cedera Fatal | Cedera Tidak Fatal | Total |
---|---|---|---|
Tidak memakai sabuk pengaman | 1,601 | 162,527 | 164,128 |
Memakai sabuk pengaman | 510 | 412,368 | 412,878 |
Total | 2,111 | 574,895 | 577,006 |
Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988
Proporsi pada No Seatbelt:
\[ \hat{p}_1 = \frac{1601}{164128} \approx 0.0097 \]
Proporsi pada Seatbelt:
\[ \hat{p}_2 = \frac{510}{412878} \approx 0.0012 \]
\[ \hat{p} = \frac{1601+510}{162527+412368} = \frac{2111}{574895} \approx 0.0036 \]
Rumus:
\[ Z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{\hat{p}(1 - \hat{p})\left(\frac{1}{n_1} + \frac{1}{n_2} \right)}} \]
\[ Z = \frac{0.0097 - 0.0012}{\sqrt{0.0036(1 - 0.0036)\left(\frac{1}{164128} + \frac{1}{412878} \right)}} \approx 48.35712 \]
x1 <- 1601
n1 <- 164128
x2 <- 510
n2 <- 412878
#Proporsi masing masing kelompok
p1 <- x1 / n1
p2 <- x2 / n2
#Proporsi gabungan
p_pool <- (x1 + x2) / (n1 + n2)
p_pool
## [1] 0.003658541
#Statistik Uji Z
Z <- (p1 - p2) / sqrt(p_pool * (1 - p_pool) * (1/n1 + 1/n2))
Z
## [1] 48.35712
# Tampilkan hasil
p1
## [1] 0.009754582
p2
## [1] 0.001235232
p_pool
## [1] 0.003658541
Z
## [1] 48.35712
Interpretasi : Jika kita menggunakan taraf signifikansi sebesar 5%, kita menerima H0 karena nilai |Z| yaitu 48.35712 > dari 1.96 yang berarti bahwa ada perbedaan proporsi antar kelompok.
Uji asosiasi dalam tabel kontingensi 2×2 bertujuan untuk mengukur hubungan antara dua variabel kategori. Untuk setiap uji asosiasi, hipotesis yang diuji adalah :
Tiga ukuran utama dalam uji asosiasi adalah:
Risk difference mengukur perbedaan absolut dalam probabilitas kejadian antara dua kelompok.
Rumus : \[ RD = \frac{a}{a + b} - \frac{c}{c + d} \] Standar Error untuk RD \[ SE(RD) = \sqrt{ \frac{p_1(1 - p_1)}{n_1} + \frac{p_2(1 - p_2)}{n_2} } \]
Statistik uji Z untuk RD \[ Z = \frac{RD}{SE(RD)} \]
Perbandingan antara risiko kejadian di kelompok 1 dibanding kelompok 2.
Rumus:
\[ RR = \frac{a / (a + b)}{c / (c + d)} \]
Standar Error untuk RR \[ SE(\log(RR)) = \sqrt{ \frac{1 - p_1}{n_1 p_1} + \frac{1 - p_2}{n_2 p_2} } \] Statistik uji Z untuk RR \[ Z = \frac{\log(RR)}{SE(\log(RR))} \]
Perbandingan peluang kejadian antara dua kelompok
Rumus: \[ OR = \frac{a \cdot d}{b \cdot c} \] \[ SE(\log(OR)) = \sqrt{ \frac{1}{a} + \frac{1}{b} + \frac{1}{c} + \frac{1}{d} } \]
Statistik uji Z untuk OR \[ Z = \frac{\log(OR)}{SE(\log(OR))} \]
Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.
Penggunaan Sabuk Pengaman | Cedera Fatal | Cedera Tidak Fatal | Total |
---|---|---|---|
Tidak memakai sabuk pengaman | 1,601 | 162,527 | 164,128 |
Memakai sabuk pengaman | 510 | 412,368 | 412,878 |
Total | 2,111 | 574,895 | 577,006 |
Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988
\[ RD = \hat{p}_1 - \hat{p}_2 = 0.00976 - 0.00124 = 0.00852 \]
\[ SE(RD) = \sqrt{ \frac{0.00976(1 - 0.00976)}{164128} + \frac{0.00124(1 - 0.00124)}{412878} } \] \[ = \sqrt{5.89 \times 10^{-8} + 3.01 \times 10^{-9}} = \sqrt{6.19 \times 10^{-8}} = 0.000249 \]
\[ Z = \frac{0.00852}{0.000249} \approx 34.22 \]
\[ RR = \frac{0.00976}{0.00124} = 7.87 \]
Standard Error log(RR)
\[ SE(\log(RR)) = \sqrt{ \left( \frac{1}{1601} - \frac{1}{164128} \right) + \left( \frac{1}{510} - \frac{1}{412878} \right) } \]
\[ = \sqrt{(0.0006246 - 0.0000061) + (0.0019608 - 0.0000024)} = \sqrt{0.0006185 + 0.0019584} \] \[ = \sqrt{0.0025769} = 0.05076 \] Statistik uji Z untuk RR
\[ Z_{RR} = \frac{\log(7.87)}{0.05076} = \frac{2.06}{0.05076} \approx 40.58 \]
\[ OR = \frac{a \cdot d}{b \cdot c} = \frac{1601 \cdot 412368}{162527 \cdot 510} = \frac{660242368}{82888770} = 7.97 \] Standard Error log(OR)
\[ SE(\log(OR)) = \sqrt{ \frac{1}{1601} + \frac{1}{162527} + \frac{1}{510} + \frac{1}{412368} } = \sqrt{0.0006246 + 0.0000062 + 0.0019608 + 0.0000024} = \sqrt{0.0025939} = 0.05093 \] Statistik uji Z untuk OR
\[ Z_{OR} = \frac{\log(7.97)}{0.05093} = \frac{2.076}{0.05093} \approx 40.77 \]
Uji independensi digunakan untuk menentukan apakah ada hubungan statistik antara dua variabel kategorikal.
Uji Chi-Square digunakan untuk menguji apakah ada hubungan antara dua variabel kategorikal.
Rumus Chi-Square
\[ \chi^2 = \sum \frac{(O - E)^2}{E} \]
di mana:
\[ E_{ij} = \frac{R_i \times C_j}{N} \] dengan:
Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.
Penggunaan Sabuk Pengaman | Cedera Fatal | Cedera Tidak Fatal | Total |
---|---|---|---|
Tidak memakai sabuk pengaman | 1,601 | 162,527 | 164,128 |
Memakai sabuk pengaman | 510 | 412,368 | 412,878 |
Total | 2,111 | 574,895 | 577,006 |
Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988
Perhitungan Manual
Langkah 1: Hitung nilai yang diharapkan (E)
Nilai yang diharapkan \(E_{ij}\) dihitung dengan rumus:
\[ E_{ij} = \frac{R_i \times C_j}{N} \]
dimana:
Menghitung nilai yang diharapkan:
Gunakan rumus:
\[ E_{ij} = \frac{\text{Total baris } i \times \text{Total kolom } j}{\text{Total keseluruhan}} \]
Untuk Tidak Pakai Sabuk, Cedera Fatal (\(E_{11}\)):
\[ E_{11} = \frac{164128 \times 2111}{577006} = \frac{346555408}{577006} \approx 600.2 \]
Untuk Tidak Pakai Sabuk, Tidak Fatal (\(E_{12}\)):
\[ E_{12} = \frac{164128 \times 574895}{577006} \approx 163527.8 \]
Untuk Pakai Sabuk, Cedera Fatal (\(E_{21}\)):
\[ E_{21} = \frac{412878 \times 2111}{577006} \approx 1510.8 \]
Untuk Pakai Sabuk, Tidak Fatal (\(E_{22}\)):
\[ E_{22} = \frac{412878 \times 574895}{577006} \approx 411367.2 \]
Hitung Chi-Square (\(\chi^2\))
Gunakan rumus:
\[ \chi^2 = \sum \frac{(O - E)^2}{E} \]
Untuk Tidak Pakai Sabuk, Fatal:
\[ \frac{(1601 - 600.2)^2}{600.2} = \frac{1000800.64}{600.2} \approx 1667.9 \]
Untuk Tidak Pakai Sabuk, Tidak Fatal:
\[ \frac{(162527 - 163527.8)^2}{163527.8} = \frac{1008010.24}{163527.8} \approx 6.2 \]
Untuk Pakai Sabuk, Fatal:
\[ \frac{(510 - 1510.8)^2}{1510.8} = \frac{1008010.24}{1510.8} \approx 667.2 \]
Untuk Pakai Sabuk, Tidak Fatal:
\[ \frac{(412368 - 411367.2)^2}{411367.2} = \frac{1000800.64}{411367.2} \approx 2.4 \]
Total Chi-Square:
\[ \chi^2 = 1667.9 + 6.2 + 667.2 + 2.4 = 2343.7 \]
Perhitungan menggunakan R
# Membuat data tabel 2x2
data <- matrix(c(1601, 162527, 510, 412368), nrow = 2, byrow = TRUE)
# Menambahkan nama baris dan kolom
rownames(data) <- c("Tidak Memakai Sabuk", "Memakai Sabuk")
colnames(data) <- c("Cedera Fatal", "Cedera Tidak Fatal")
# Menampilkan data tabel
data
## Cedera Fatal Cedera Tidak Fatal
## Tidak Memakai Sabuk 1601 162527
## Memakai Sabuk 510 412368
# Menghitung uji chi-square
chi_square_test <- chisq.test(data)
# Menampilkan hasil uji chi-square
chi_square_test
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: data
## X-squared = 2336.1, df = 1, p-value < 2.2e-16
Interpretasi: Dengan taraf signifikansi 5%, didapatkan p-value < 0.05, maka kita menolak H0 (hipotesis nol). Artinya, ada hubungan yang signifikan secara statistik antara penggunaan sabuk pengaman dan tingkat fatalitas kecelakaan
Partisi Chi-Square Partisi Chi-Square digunakan untuk mengidentifikasi kategori mana dalam tabel kontingensi yang bertanggung jawab atas hubungan yang signifikan. Jika uji Chi-Square pada tabel kontingensi I×Jsignifikan, maka partisi Chi-Square memungkinkan kita untuk menguraikan efek hubungan dalam subkelompok yang lebih kecil.
Contoh Kasus:
Identifikasi Partai Berdasarkan Ras (GSS 1991)
Ras | Demokrat | Independen | Republik |
---|---|---|---|
Black | 103 | 15 | 11 |
White | 341 | 105 | 405 |
Total | 444 | 120 | 416 |
Sumber: 1991 General Social Survey, National Opinion Research Center
Perhitungan Manual Langkah 1: Menghitung Frekuensi yang Diharapkan
Frekuensi yang diharapkan dihitung menggunakan rumus:
\[ E_{ij} = \frac{(Total \, \text{Baris}) \times (Total \, \text{Kolom})}{Grand \, \text{Total}} \]
\[ E_{1,1} = \frac{129 \times 444}{980} = 58.48 \]
\[ E_{1,2} = \frac{129 \times 120}{980} = 15.79 \]
\[ E_{1,3} = \frac{129 \times 416}{980} = 54.73 \]
\[ E_{2,1} = \frac{851 \times 444}{980} = 385.52 \]
\[ E_{2,2} = \frac{851 \times 120}{980} = 104.20 \]
\[ E_{2,3} = \frac{851 \times 416}{980} = 361.27 \]
Langkah 2: Menghitung Statistik Chi-Square
Sekarang kita dapat menghitung statistik chi-square \[\chi^2\] dengan rumus:
\[ \chi^2 = \sum \frac{(O_{ij} - E_{ij})^2}{E_{ij}} \]
Dimana \(O_{ij}\) adalah frekuensi yang diamati dan \(E_{ij}\) adalah frekuensi yang diharapkan.
Dimana \(O_{ij}\) adalah frekuensi yang diamati dan \(E_{ij}\) adalah frekuensi yang diharapkan. Berdasarkan data yang ada, kita menghitung chi-square untuk setiap sel:
\[ \chi^2 = \frac{(103 - 58.48)^2}{58.48} + \frac{(15 - 15.79)^2}{15.79} + \frac{(11 - 54.73)^2}{54.73}+ \] \[ \frac{(341 - 385.52)^2}{385.52} + \frac{(105 - 104.20)^2}{104.20} + \frac{(405 - 361.27)^2}{361.27} \]
\[ = 33.89 + 0.0395 + 34.95 + 5.14 + 0.0061 + 5.29 = 79.32 \]
\[ \chi^2 \approx 79.32 \]
\[ \chi^2 = 33.89 + 0.0395 + 34.95 + 5.14 + 0.0061 + 5.29 = 79.32 \]
Jadi, berdasarkan perhitungan manual:
\[ \chi^2 \approx 79.32 \]
Perhitungan Dengan R:
# Membuat matriks data
data_matrix <- matrix(c(103, 15, 11,
341, 105, 405),
nrow = 2, byrow = TRUE)
colnames(data_matrix) <- c("Demokrat", "Independen", "Republik")
rownames(data_matrix) <- c("Black", "White")
# Melakukan uji chi-square
chi_square_test <- chisq.test(data_matrix)
# Menampilkan hasil uji
chi_square_test
##
## Pearson's Chi-squared test
##
## data: data_matrix
## X-squared = 79.431, df = 2, p-value < 2.2e-16
Interpretasi: Karena nilai \(\chi^2\) sangat besar dan jauh melebihi nilai kritis untuk df = 2 pada taraf signifikansi 0.05, maka kita menolak \(H_0\). Ini menunjukkan bahwa ada hubungan signifikan antara ras dan identifikasi partai.
Uji Likelihood Ratio (G²) adalah alternatif dari uji chi-square yang digunakan untuk menguji hipotesis independensi dalam tabel kontingensi \(I \times J\). Statistik uji ini diberikan oleh:
\[ G^2 = 2 \sum_{i} \sum_{j} n_{ij} \ln \left( \frac{n_{ij}}{\hat{\mu}_{ij}} \right) \]
Dimana:
- \(n_{ij}\) adalah frekuensi observasi dalam tabel kontingensi.
\(\hat{\mu}_{ij} = n \cdot p_i \cdot p_j\) adalah frekuensi yang diharapkan.
\(p_i\) adalah proporsi untuk baris \(i\).
- \(p_j\) adalah proporsi untuk kolom \(j\).
Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.
Penggunaan Sabuk Pengaman | Cedera Fatal | Cedera Tidak Fatal | Total |
---|---|---|---|
Tidak memakai sabuk pengaman | 1,601 | 162,527 | 164,128 |
Memakai sabuk pengaman | 510 | 412,368 | 412,878 |
Total | 2,111 | 574,895 | 577,006 |
Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988
Langkah 1: Hitung Frekuensi Ekspektasi
Gunakan rumus: \[ \hat{\mu}_{ij} = \frac{(\text{Total Baris}_i) \times (\text{Total Kolom}_j)}{\text{Grand Total}} \]
\[ \begin{aligned} \hat{\mu}_{11} &= \frac{164128 \times 2111}{577006} = 599.78 \\ \hat{\mu}_{12} &= \frac{164128 \times 574895}{577006} = 163528.22 \\ \hat{\mu}_{21} &= \frac{412878 \times 2111}{577006} = 1511.22 \\ \hat{\mu}_{22} &= \frac{412878 \times 574895}{577006} = 411366.78 \\ \end{aligned} \] Langkah 2: Hitung Statistik Uji G²
\[ \begin{aligned} G^2 &= 2 \sum O_{ij} \ln \left( \frac{O_{ij}}{E_{ij}} \right) \\ &= 2 \bigg[ 1601 \ln \left( \frac{1601}{599.78} \right) + 162527 \ln \left( \frac{162527}{163528.22} \right) + \\ &\quad 510 \ln \left( \frac{510}{1511.22} \right) + 412368 \ln \left( \frac{412368}{411366.78} \right) \bigg] \end{aligned} \]
Setelah dihitung:
\[ \begin{aligned} G^2 &= 2 \times (1693.26 + (-6.14) + (-1037.45) + 2.43) \\ &= 2 \times 652.10 = 1304.20 \end{aligned} \]
Langkah 3: Bandingkan dengan Distribusi Chi-Square
Derajat bebas:
\[ (I - 1)(J - 1) = (2 - 1)(2 - 1) = 1 \]
Nilai kritis \(\chi^2\) untuk \(\alpha = 0.05\) dan df = 1 adalah 3.841.
Perhitungan Dengan R
# Data Observasi
observed <- matrix(c(1601, 162527,
510, 412368),
nrow = 2, byrow = TRUE)
rownames(observed) <- c("Tidak Pakai Sabuk", "Pakai Sabuk")
colnames(observed) <- c("Cedera Fatal", "Tidak Fatal")
# Total per baris dan kolom
row_totals <- rowSums(observed)
col_totals <- colSums(observed)
grand_total <- sum(observed)
# Frekuensi harapan
expected <- outer(row_totals, col_totals) / grand_total
# Hitung G²
G2 <- 2 * sum(observed * log(observed / expected), na.rm = TRUE)
# Derajat bebas
df <- (nrow(observed) - 1) * (ncol(observed) - 1)
# Nilai kritis
chi_critical <- qchisq(0.95, df)
# Output hasil
cat("Pearson's Chi Square =", round(G2, 2), "\n")
## Pearson's Chi Square = 2041.16
cat("Derajat bebas =", df, "\n")
## Derajat bebas = 1
cat("Nilai kritis chi-square (0.05) =", round(chi_critical, 3), "\n")
## Nilai kritis chi-square (0.05) = 3.841
if (G2 > chi_critical) {
cat("Kesimpulan: Tolak H0 - Ada hubungan signifikan.\n")
} else {
cat("Kesimpulan: Gagal tolak H0 - Tidak ada hubungan signifikan.\n")
}
## Kesimpulan: Tolak H0 - Ada hubungan signifikan.
Interpretasi: Karena \(G^2 = 1304.20 > 3.841\), maka kita menolak hipotesis nol, artinya terdapat bukti yang cukup bahwa penggunaan sabuk pengaman berpengaruh signifikan terhadap tingkat fatalitas kecelakaan.
Tabel kontingensi tiga arah adalah perpanjangan dari tabel kontingensi dua arah yang digunakan untuk menganalisis hubungan antara tiga variabel kategori secara simultan. Dalam banyak situasi, hubungan antara dua variabel (misalnya X dan Y) dapat dipengaruhi oleh variabel ketiga Z, yang disebut sebagai variabel kontrol. Tabel ini sering digunakan dalam analisis data ketika ada potensi faktor pengganggu yang dapat mempengaruhi hubungan antara dua variabel utama. Misalnya:
Tabel kontingensi tiga arah dapat dibagi menjadi dua jenis:
Tabel Parsial: Tabel yang menyajikan hubungan antara X dan Y pada setiap kategori Z. Ini memungkinkan analisis hubungan bersyarat, dimana efek Z dikendalikan.
Tabel Marginal: Tabel yang diperoleh dengan mengabaikan Z, yaitu dengan menjumlahkan semua kategori Z. Ini memberikan gambaran umum hubungan antara X dan Y tanpa mempertimbangkan Z, yang bisa menyebabkan distorsi interpretasi, seperti dalam Simpson’s Paradox
Kegunaan tabel marginal adalah untuk melihat pola asosiasi secara agregat, tetapi sering kali mengabaikan efek kovariat yang dapat memberikan pemahaman lebih mendalam. Oleh karena itu, analisis yang mempertimbangkan tabel parsial biasanya lebih ditentukan untuk mendapatkan kesimpulan yang lebih akurat dalam penelitian ilmiah.
Tabel parsial adalah tabel yang mengelompokkan𝑋dan𝑌berdasarkan setiap level 𝑍, sedangkan tabel marginal adalah tabel yang mengabaikan𝑍,dengan menjumlahkan data dari semua level𝑍.
Tabel berikut menunjukkan data mengenai efektivitas obat AZT (zidovudine) dalam memperlambat perkembangan gejala AIDS. Studi ini melibatkan 338 veteran yang rasnya dibagi menjadi dua yaitu ras kulit putih dan ras kulit hitam yang telah terinfeksi virus AIDS dan mulai menunjukkan tanda-tanda penurunan sistem kekebalan tubuh.
Race | AZT Use | Yes (Symptoms) | No (Symptoms) |
---|---|---|---|
White | Yes | 14 | 93 |
White | No | 32 | 81 |
Black | Yes | 11 | 52 |
Black | No | 12 | 43 |
Sumber: New York Times, Feb 15, 1991.
Tabel Frekuensi Parsial untuk Z (Race) = White
AZT Use | Yes (Symptoms) | No (Symptoms) |
---|---|---|
Yes | 14 | 93 |
No | 32 | 81 |
Tabel Frekuensi Parsial untuk Z (Race) = White menunjukkan hubungan antara penggunaan AZT dan perkembangan gejala AIDS pada kelompok ras kulit putih.
Tabel Frekuensi Parsial untuk Z (Race) = Black
AZT Use | Yes (Symptoms) | No (Symptoms) |
---|---|---|
Yes | 11 | 52 |
No | 12 | 43 |
Tabel Frekuensi Parsial untuk Z (Race) = Black menunjukkan hubungan antara penggunaan AZT dan perkembangan gejala AIDS pada kelompok ras kulit hitam.
Dengan R:
# Membuat data dengan array
data_azt <- array(
c(14, 93, 32, 81, 11, 52, 12, 43),
dim = c(2, 2, 2),
dimnames = list(
AZT_Use = c("Yes", "No"),
Symptoms = c("Yes", "No"),
Race = c("White", "Black")
)
)
# Tampilkan array untuk memastikan datanya
data_azt
## , , Race = White
##
## Symptoms
## AZT_Use Yes No
## Yes 14 32
## No 93 81
##
## , , Race = Black
##
## Symptoms
## AZT_Use Yes No
## Yes 11 12
## No 52 43
# Ekstrak tabel frekuensi parsial berdasarkan ras
freq_parsial_white <- data_azt[, , "White"]
freq_parsial_black <- data_azt[, , "Black"]
# Tampilkan hasil frekuensi parsial
freq_parsial_white
## Symptoms
## AZT_Use Yes No
## Yes 14 32
## No 93 81
freq_parsial_black
## Symptoms
## AZT_Use Yes No
## Yes 11 12
## No 52 43
Tabel frekuensi marginal menampilkan jumlah total observasi untuk setiap variabel dengan mengabaikan variabel lainnya dalam tabel kontingensi tiga arah. Tabel ini membantu dalam memahami distribusi kategori secara agregat tanpa mempertimbangkan hubungan antarvariabel. Tabel frekuensi marginal dihitung dengan menjumlahkan frekuensi dari tabel kontingensi tiga arah berdasarkan variabel yang tersisa.
Tabel frekuensi marginal untuk Race dan Symptoms
Race | Yes | No | Total |
---|---|---|---|
White | 46 | 174 | 220 |
Black | 23 | 95 | 118 |
Total | 69 | 269 | 338 |
Tabel frekuensi marginal untuk AZT Use dan Symptoms
AZT Use | Yes | No | Total |
---|---|---|---|
Yes | 25 | 145 | 170 |
No | 44 | 124 | 168 |
Total | 69 | 269 | 338 |
Dengan R:
# Buat data array 3D: AZT Use x Symptoms x Race
data_azt <- array(
c(
14, 32, # White - AZT Yes
93, 81, # White - AZT No
11, 12, # Black - AZT Yes
52, 43 # Black - AZT No
),
dim = c(2, 2, 2),
dimnames = list(
AZT_Use = c("Yes", "No"),
Symptoms = c("Yes", "No"),
Race = c("White", "Black")
)
)
# Tampilkan array
data_azt
## , , Race = White
##
## Symptoms
## AZT_Use Yes No
## Yes 14 93
## No 32 81
##
## , , Race = Black
##
## Symptoms
## AZT_Use Yes No
## Yes 11 52
## No 12 43
# Frekuensi marginal untuk AZT Use (menjumlah semua Symptoms & Race)
freq_marginal_azt <- apply(data_azt, 1, sum)
# Frekuensi marginal untuk Symptoms (menjumlah semua AZT Use & Race)
freq_marginal_symptoms <- apply(data_azt, 2, sum)
# Tampilkan hasil
freq_marginal_azt
## Yes No
## 170 168
freq_marginal_symptoms
## Yes No
## 69 269
Tabel berikut menunjukkan data mengenai efektivitas obat AZT (zidovudine) dalam memperlambat perkembangan gejala AIDS. Studi ini melibatkan 338 veteran yang rasnya dibagi menjadi dua yaitu ras kulit putih dan ras kulit hitam yang telah terinfeksi virus AIDS dan mulai menunjukkan tanda-tanda penurunan sistem kekebalan tubuh.
Race | AZT Use | Yes (Symptoms) | No (Symptoms) |
---|---|---|---|
White | Yes | 14 | 93 |
White | No | 32 | 81 |
Black | Yes | 11 | 52 |
Black | No | 12 | 43 |
Peluang bersama didefinisikan sebagai:
\[ P(Z, X, Y) = \frac{f(Z, X, Y)}{N} \]
Contoh Perhitungan dengan Data diatas
# Load library yang diperlukan
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
# Data dalam format array 3D
data_azt <- array(
c(
14, 32, # White - AZT Yes
93, 81, # White - AZT No
11, 12, # Black - AZT Yes
52, 43 # Black - AZT No
),
dim = c(2, 2, 2),
dimnames = list(
AZT_Use = c("Yes", "No"),
Symptoms = c("Yes", "No"),
Race = c("White", "Black")
)
)
# Menampilkan tabel array 3D
print(data_azt)
## , , Race = White
##
## Symptoms
## AZT_Use Yes No
## Yes 14 93
## No 32 81
##
## , , Race = Black
##
## Symptoms
## AZT_Use Yes No
## Yes 11 52
## No 12 43
total <- sum(data_azt)
joint_prob <- data_azt/ total
ftable(joint_prob)
## Race White Black
## AZT_Use Symptoms
## Yes Yes 0.04142012 0.03254438
## No 0.27514793 0.15384615
## No Yes 0.09467456 0.03550296
## No 0.23964497 0.12721893
Interpretasi: Secara umum, data menunjukkan bahwa AZT mungkin efektif dalam menurunkan gejala AIDS, tapi belum bisa disimpulkan secara pasti tanpa analisis statistik lanjutan.
Peluang marginal didefinisikan sebagai berikut
untuk Ras: \[ P(Z) = \frac{n(Z)}{N} \]
marginal_Z <- apply(joint_prob, 3, sum)
marginal_Z
## White Black
## 0.6508876 0.3491124
Interpretasi : Proporsi Ras Putih lebih banyak dibanding ras Hitam.
Untuk AZT use: \[ P(X) = \frac{n(X)}{N} \]
marginal_X <- apply(joint_prob, 1, sum)
marginal_X
## Yes No
## 0.5029586 0.4970414
Intepretasi : Proporsi penggunaan AZT lebih banyak dibandingkan yang tidak memakai AZT
Contoh: Mencari peluang seseorang ber-ras tertentu jika ia mendapat AZT dan mengalami gejala AIDS:
\[ P(Race| AZT = Yes, Symptoms = Yes) = \frac{P(Race, AZT = Yes, Symptoms = Yes)}{P(AZT = Yes, Symptoms = Yes)} \]
# Hitung peluang bersyarat P(Race | AZT Use, Symptoms)
p_Race_given_AZT_Symptoms <- prop.table(data_azt, margin = c(1, 2))
p_Race_given_AZT_Symptoms
## , , Race = White
##
## Symptoms
## AZT_Use Yes No
## Yes 0.5600000 0.6413793
## No 0.7272727 0.6532258
##
## , , Race = Black
##
## Symptoms
## AZT_Use Yes No
## Yes 0.4400000 0.3586207
## No 0.2727273 0.3467742
Interpretasi:
Peluang seseorang ber-ras putih jika ia mendapat AZT dan mengalami gejala AIDS adalah 0.56.
Peluang seseorang ber-ras hitam jika ia mendapat AZT dan mengalami gejala AIDS adalah 0.44.
\[ RD = \frac{a}{a + b} - \frac{c}{c + d} \]
#Hitung Risk Difference (RD) untuk masing-masing ras
# White
rd_white <- (data_azt["Yes", "Yes", "White"] / sum(data_azt["Yes", , "White"])) -
(data_azt["No", "Yes", "White"] / sum(data_azt["No", , "White"]))
# Black
rd_black <- (data_azt["Yes", "Yes", "Black"] / sum(data_azt["Yes", , "Black"])) -
(data_azt["No", "Yes", "Black"] / sum(data_azt["No", , "Black"]))
# Tampilkan hasil
rd_white
## [1] -0.1523447
rd_black
## [1] -0.04357864
Interpretasi: AZT mengurangi risiko berkembangnya gejala AIDS sebesar 15.2% pada ras White. AZT mengurangi risiko sebesar 4.36% pada ras Black, tapi efeknya lebih kecil dibanding pada ras White.
\[ RR = \frac{a / (a + b)}{c / (c + d)} \]
#Hitung Risk Relative (RR) untuk masing-masing ras
# White
rr_white <- (data_azt["Yes", "Yes", "White"] / sum(data_azt["Yes", , "White"])) /
(data_azt["No", "Yes", "White"] / sum(data_azt["No", , "White"]))
# Black
rr_black <- (data_azt["Yes", "Yes", "Black"] / sum(data_azt["Yes", , "Black"])) /
(data_azt["No", "Yes", "Black"] / sum(data_azt["No", , "Black"]))
# Tampilkan hasil
rr_white
## [1] 0.4620327
rr_black
## [1] 0.8002646
Interpretasi: Pada kelompok ras White, risiko berkembangnya gejala AIDS bagi yang menggunakan AZT adalah sekitar 46% dari risiko pada yang tidak menggunakan AZT. Pada kelompok ras Black, risiko berkembangnya gejala AIDS bagi yang menggunakan AZT adalah sekitar 80% dari risiko pada yang tidak menggunakan AZT.
\[ OR = \frac{a \cdot d}{b \cdot c} \]
# Menghitung Odds Ratio untuk masing masing ras
or_white <- (data_azt["Yes", "Yes", "White"] / data_azt["Yes", "No", "White"]) /
(data_azt["No", "Yes", "White"] / data_azt["No", "No", "White"])
or_black <- (data_azt["Yes", "Yes", "Black"] / data_azt["Yes", "No", "Black"]) /
(data_azt["No", "Yes", "Black"] / data_azt["No", "No", "Black"])
or_white
## [1] 0.3810484
or_black
## [1] 0.7580128
Interpretasi: Di antara individu dengan ras putih, peluang (odds) mengalami gejala AIDS pada yang menggunakan AZT adalah sekitar 38% dari odds mengalami gejala pada yang tidak menggunakan AZT Di antara individu dengan ras hitam, peluang (odds) mengalami gejala AIDS pada yang menggunakan AZT adalah sekitar 76% dari odds mengalami gejala pada yang tidak menggunakan AZT
Tabel kontingensi tiga arah digunakan untuk menganalisis hubungan antara dua variabel kategorik dengan mempertimbangkan variabel kontrol. Contohnya adalah hubungan antara Penggunaan AZT (X) dan perkembangan gejala AIDS (Y), dengan variabel kontrol ras (Z)
Tabel ini terdiri dari beberapa tabel parsial (2×2) untuk setiap tingkat Z, serta tabel marginal yang mengabaikan Z. Ukuran asosiasi yang digunakan adalah odds ratio.
Jika odds ratio parsial relatif konstan, kita dapat menghitung odds ratio bersama menggunakan estimasi Mantel-Haenszel.
Independensi bersyarat adalah konsep penting dalam analisis tabel kontingensi tiga arah. Ini merujuk pada kondisi di mana dua variabel, 𝑋 dan 𝑌, independen dalam setiap level variabel ketiga, 𝑍. Pengujian independensi bersyarat dilakukan dengan metode statistik seperti uji Cochran-Mantel-Haenszel (CMH).
Definisi Independensi Bersyarat
Independensi Bersyarat: Dua variabel, 𝑋 dan 𝑌, dikatakan independen bersyarat terhadap variabel ketiga, 𝑍, jika rasio odds mereka dalam setiap strata 𝑍 sama dengan 1.
Metode Cochran-Mantel-Haenszel Tujuan Uji CMH Uji Cochran-Mantel-Haenszel (CMH) digunakan untuk menguji hubungan antara dua variabel kategori dengan mempertimbangkan efek dari variabel perancu (confounder).
Statistik uji Cochran-Mantel-Haenszel (CMH) dirumuskan sebagai:
\[ CMH = \frac{\sum_k \left( n_{1ik} - \mu_{1ik} \right)^2}{\sum_k \text{var}(n_{1ik})} \]
Keterangan
\[ \mu_{1ik} = E(n_{1ik}) = \frac{n_{1.k} \cdot n_{1ik}}{n_{.k}} \]
Varians dari \(n_{1ik}\) diberikan oleh:
\[ \text{var}(n_{1ik}) = \frac{n_{1.k} \cdot n_{2.k} \cdot n_{1ik} \cdot n_{2.k}}{n_{2.k}^2 \cdot (n_{.k} - 1)} \] Contoh Perhitungan:
Race: White
Symptoms (Yes) | Symptoms (No) | Total | |
---|---|---|---|
AZT Yes | 14 | 93 | 107 |
AZT No | 32 | 81 | 113 |
Race: Black
Symptoms (Yes) | Symptoms (No) | Total | |
---|---|---|---|
AZT Yes | 11 | 52 | 63 |
AZT No | 12 | 43 | 55 |
Perhitungan Manual:
\[ \mu_{1ik} = \frac{n_{i+k} \cdot n_{+ik}}{n_{++k}} \]
Untuk Ras Putih (k=1):
\[ \mu_{111} = \frac{107 \times 46}{220} = 22.38 \]
Untuk Ras Hitam (k=2):
\[ \mu_{112} = \frac{63 \times 23}{118} = 12.28 \]
\[ \text{Var}(n_{1ik}) = \frac{n_{i+k} \cdot n_{0+k} \cdot n_{+ik} \cdot n_{+0k}}{n_{++k}^2 \cdot (n_{++k} - 1)} \]
Untuk Ras Putih (k=1):
\[ \text{Var}(n_{111}) = \frac{107 \cdot 113 \cdot 46 \cdot 174}{220^2 \cdot 219} \approx 9.41 \]
Untuk Ras Hitam (k=2):
\[ \text{Var}(n_{112}) = \frac{63 \cdot 55 \cdot 23 \cdot 95}{118^2 \cdot 117} \approx 5.25 \]
\[ X^2_{CMH} = \frac{(n_{111} - \mu_{111})^2 + (n_{112} - \mu_{112})^2}{\text{Var}(n_{111}) + \text{Var}(n_{112})} \]
\[ X^2_{CMH} = \frac{(14 - 22.38)^2 + (11 - 12.28)^2}{9.41 + 5.25} = \frac{74.16}{14.66} \approx 5.06 \]
Interpretasi: Dengan taraf signifikansi 5%, nilai kritis \(\chi^2\) adalah 3.841. Karena \(X^2_{CMH} = 5.06 > 3.841\), maka hipotesis nol ditolak. Artinya, ada hubungan signifikan antara penggunaan AZT dan perkembangan gejala AIDS setelah mengontrol variabel ras.
Perhitungan dengan R:
# Data dalam format array 3D
data_azt <- array(
c(
14, 32, # White - AZT Yes
93, 81, # White - AZT No
11, 12, # Black - AZT Yes
52, 43 # Black - AZT No
),
dim = c(2, 2, 2),
dimnames = list(
AZT_Use = c("Yes", "No"),
Symptoms = c("Yes", "No"),
Race = c("White", "Black")
)
)
# Menampilkan tabel
print(data_azt)
## , , Race = White
##
## Symptoms
## AZT_Use Yes No
## Yes 14 93
## No 32 81
##
## , , Race = Black
##
## Symptoms
## AZT_Use Yes No
## Yes 11 52
## No 12 43
cmh_base <- mantelhaen.test(data_azt, correct = FALSE)
cmh_base
##
## Mantel-Haenszel chi-squared test without continuity correction
##
## data: data_azt
## Mantel-Haenszel X-squared = 6.7624, df = 1, p-value = 0.00931
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 0.2818621 0.8414208
## sample estimates:
## common odds ratio
## 0.4869955
Interpretasi : Dikarenakan p-value < 0.05, maka keputusan uji adalah menolak H₀, yang berarti terdapat hubungan yang signifikan antara penggunaan AZT dengan perkembangan gejala AIDS setelah mengontrol variabel ras. Dengan kata lain, AZT secara signifikan berhubungan dengan timbul atau tidaknya gejala AIDS, terlepas dari perbedaan ras peserta dalam studi ini
Generalized Linear Model (GLM) merupakan perluasan dari model regresi linear klasik. GLM memungkinkan digunakan untuk memodelkan data dimana variabel respons tidak berdistribusi normal dan atau hubungan antara variabel prediktor dengan rata-rata dari variabel respons tidak linear. GLM terdiri dari tiga komponen utama:
Distribusi termasuk dalam exponential family jika dapat ditulis dalam bentuk:
\[ f(y; \theta, \phi) = \exp \left\{ \frac{y\theta - b(\theta)}{\phi} + c(y, \phi) \right\} \]
Contoh Distribusi yang Termasuk Exponential Family:
Contoh Pembuktian: Distribusi Binomial
Fungsi probabilitas distribusi binomial adalah:
\[ P(Y = y) = \binom{n}{y} \pi^y (1 - \pi)^{n - y} \]
Tuliskan ulang dalam bentuk exponential family:
\[ P(Y = y) = \exp \left\{ \log \binom{n}{y} + y \log \left( \frac{\pi}{1 - \pi} \right) + n \log (1 - \pi) \right\} \]
Dengan substitusi:
Kesimpulan : Dengan bentuk tersebut, maka distribusi binomial termasuk dalam keluarga exponential family.
Regresi Logistik memiliki kemiripan dengan regresi linear, yaitu menggabungkan nilai-nilai input secara linear dengan bobot (koefisien) untuk menghasilkan sebuah prediksi. Bedanya, regresi logistik membatasi hasil prediksinya agar berada dalam bentuk nilai biner 0 atau 1 dengan bantuan fungsi sigmoid sebagai fungsi aktivasi. Fungsi ini mengubah output menjadi nilai antara 0 sampai 1 dan membentuk kurva seperti huruf S.
Model ini digunakan untuk menganalisis hubungan antara satu atau lebih variabel independen dan mengklasifikasikan data ke dalam kelas-kelas yang bersifat diskrit. Artinya, regresi logistik cocok digunakan untuk masalah klasifikasi, terutama yang hanya memiliki dua kemungkinan kategori (dikenal dengan klasifikasi biner).
Contohnya, angka 0 bisa merepresentasikan gagal, sedangkan angka 1 merepresentasikan sukses.
Beberapa contoh penerapan regresi logistik dalam klasifikasi biner antara lain:
-Memprediksi kelulusan siswa: Model regresi logistik dapat digunakan untuk memperkirakan apakah seorang siswa akan lulus atau tidak berdasarkan variabel-variabel seperti kehadiran, nilai tugas, partisipasi kelas, dan nilai ujian.
-Deteksi penipuan transaksi keuangan: Dalam dunia perbankan dan keuangan, regresi logistik dapat membantu mendeteksi apakah sebuah transaksi tergolong normal atau mencurigakan (fraud) dengan menganalisis pola transaksi seperti frekuensi, jumlah uang, lokasi transaksi, dan waktu.
-Memprediksi hasil diagnosis penyakit: Dalam dunia medis, model ini bisa digunakan untuk menentukan apakah pasien menderita suatu penyakit (misalnya diabetes atau kanker) berdasarkan variabel seperti usia, tekanan darah, kadar gula, dan riwayat keluarga.
Keunggulan Utama Regresi Logistik yaitu mudah diterapkan dalam Machine Learning. Regresi logistik sangat cocok digunakan dalam Machine Learning karena proses pelatihan (training) dan pengujian (testing) nya cukup sederhana. Model ini belajar dari pola-pola dalam data input dan menghubungkannya dengan output (label). Karena tidak membutuhkan komputasi tinggi, regresi logistik tergolong mudah untuk diimplementasikan, dipahami, dan dilatih dibanding metode lain dalam Machine Learning.
Cocok untuk data yang bisa dipisahkan secara linear Jika dua kelas data bisa dipisahkan dengan garis lurus (linear), maka regresi logistik akan sangat efektif dalam mengklasifikasikan data tersebut ke dalam dua kelompok berbeda. Dalam konteks ini, variabel target (y) hanya memiliki dua nilai, sehingga model bisa menentukan batas yang jelas di antara keduanya.
Memberikan wawasan yang bermakna Regresi logistik juga bisa menunjukkan seberapa besar pengaruh suatu variabel terhadap hasil akhir, melalui koefisiennya. Selain itu, model ini juga menunjukkan apakah hubungan antar variabel bersifat positif atau negatif.
Persamaan dan Asumsi dalam Regresi Logistik Model ini menggunakan fungsi sigmoid, yaitu fungsi matematika berbentuk kurva S yang mengubah nilai output ke dalam rentang antara 0 dan 1. Jika hasil dari fungsi sigmoid (yang merepresentasikan probabilitas) lebih tinggi dari ambang batas tertentu, maka model akan mengklasifikasikan data sebagai bagian dari suatu kelas. Sebaliknya, jika di bawah ambang tersebut, maka data dianggap tidak termasuk dalam kelas itu.
Selain itu, jika keluaran dari fungsi sigmoid (yaitu probabilitas yang diperkirakan) lebih besar dari ambang batas yang telah ditentukan dalam grafik, maka model akan memprediksi bahwa suatu observasi termasuk dalam kelas tertentu. Sebaliknya, jika nilai probabilitas tersebut lebih kecil dari ambang batas, model akan memprediksi bahwa observasi tersebut tidak termasuk ke dalam kelas tersebut
Sebagai contoh:
-Jika hasil fungsi sigmoid lebih dari 0,5, maka output dianggap sebagai 1 (kelas positif). -Jika hasilnya kurang dari 0,5, maka output diklasifikasikan sebagai 0 (kelas negatif). -Jika grafik menuju ke arah negatif secara ekstrem, maka nilai prediksi y akan menjadi 0,dan sebaliknya.
Fungsi sigmoid digunakan dalam regresi logistik untuk mengubah nilai prediksi menjadi probabilitas. Rumus dari fungsi sigmoid yaitu: \[ f(x) = \frac{1}{1 + e^{-x}} \]
# Simulasi data untuk regresi logistik
set.seed(42)
n <- 100
x <- seq(-4, 4, length.out = n)
log_odds <- -0.5 + 1.5 * x
prob <- 1 / (1 + exp(-log_odds))
y <- rbinom(n, 1, prob)
# Buat data frame
data <- data.frame(x=x, y=y, prob = prob)
Plot Kurva Sigmoid
# Visualisasi menggunakan base R plot
# Visualisasi menggunakan base R
# Plot the data
plot(x, y,
pch = 16,
col = "gray60",
xlab = "X",
ylab = "Y / Probabilitas",
main = "Simulasi Regresi Logistik dengan Kurva Sigmoid")
lines(x, prob, col = "blue", lwd = 2)
abline(h = 0.5, col = "red", lty = 2)
legend("topleft", legend = c("Data Biner (0/1)", "Kurva Logistik", "Ambang 0.5"),
col = c("gray60", "blue", "red"),
pch = c(16, NA, NA),
lty = c(NA, 1, 2),
lwd = c(NA, 2, 1),
pt.cex = 1.5, bty = "n")
Kurva sigmoid dalam regresi logistik menunjukkan hubungan non-linear antaravariabel prediktor dan probabilitas output. Pendekatan ini efektif untuk klasifikasi biner seperti deteksi penyakit, kelulusan siswa, dan prediksi ya/tidak.
Fungsi Link (Logit): Fungsi link logit dapat dinyatakan dalam bentuk berikut:
\[ g(\mu) = \log \left( \frac{\mu}{1 - \mu} \right) \]
Model Regresi Logistik:
\[ \log \left( \frac{\mu}{1 - \mu} \right) = X\beta \]
Fungsi Inverse Link:
\[ \mu = \frac{\exp(X\beta)}{1 + \exp(X\beta)} \]
Estimasi Parameter
Metode estimasi parameter pada GLM umumnya menggunakan Maximum Likelihood Estimation (MLE).
Log-likelihood fungsi untuk regresi logistik:
\[ l(\beta) = \sum_{i=1}^{n} \left[ y_i \log(\pi_i) + (1 - y_i) \log(1 - \pi_i) \right] \] Dengan:
\[ \pi_i = \frac{\exp(X\beta)}{1 + \exp(X\beta)} \]
Contoh Kasus dengan R
Misalkan kita memiliki:
#Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)
#Buat data biner
game <- rep(1:23, Attempts)
result <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))
#Gabung jadi data frame
data <- data.frame(Game = game, Result = result)
#Definisikan x dan y
x <- data$Game
y <- data$Result
head(data)
Estimasi Regresi Logistik Estimasi parameter model regresi logistik dapat menggunakan ‘glm’ function
model <- glm(y ~ x, data = data, family = binomial)
summary(model)
##
## Call:
## glm(formula = y ~ x, family = binomial, data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.02085 0.26206 0.080 0.937
## x -0.01546 0.01845 -0.838 0.402
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 408.06 on 295 degrees of freedom
## Residual deviance: 407.35 on 294 degrees of freedom
## AIC: 411.35
##
## Number of Fisher Scoring iterations: 3
Interpretasi Koefisien
• Intercept ((Intercept)): nilai log-odds saat x=0
• Koefisien x: perubahan log-odds untuk setiap satu satuan peningkatan pada x
• Nilai p: menunjukkan signifikansi statistik dari prediktor
Koefisien log-odds dapat ditransformasikan ke odds ratio:
exp(coef(model))
## (Intercept) x
## 1.0210698 0.9846578
Prediksi dan Visualisasi Kurva Logit
library(ggplot2)
data$pred <- predict(model, type = "response")
ggplot(data, aes(x=x, y=y)) +
geom_point(alpha = 0.5, color = "gray40") +
geom_line(aes(y=pred), color = "blue", linewidth = 1.5) +
labs(title = "Kurva Logit pada Regresi Logistik",
x="X (Prediktor)",
y="Probabilitas / Respons") +
theme_minimal()
Evaluasi Model
Kita dapat melihat akurasi model berdasarkan klasifikasi threshold 0.5
# Klasifikasi dan confusion matrix
data$pred_class <- ifelse(data$pred > 0.5, 1, 0)
table(Predicted = data$pred_class, Actual = y)
## Actual
## Predicted 0 1
## 0 160 131
## 1 1 4
GLM adalah kerangka model fleksibel untuk berbagai jenis data dan distribusi. Regresi logistik merupakan salah satu contoh penting dari GLM, sangat berguna dalam analisis data kategorik biner. Estimasi parameter dilakukan melalui metode MLE dan dapat diselesaikan secara efisien dengan fungsi glm di R.
Regresi Poisson digunakan ketika variabel respons adalah data cacah (count data),yaitu bilangan bulat non negatif. Model ini merupakan bagian dari Generalized Linear Model(GLM) dengan asumsi bahwa distribusi variabel respons adalah distribusi Poisson.
Distribusi Poisson memiliki fungsi probabilitas:
\[ P(Y = y) = \frac{e^{-\lambda} \lambda^y}{y!} \]
Distribusi Poisson dapat dituliskan dalam bentuk exponential family sebagai berikut:
\[ f(y; \theta) = \exp \left\{ y \log(\lambda) - \lambda - \log(y!) \right\} \]
dengan:
Maka distribusi Poisson termasuk dalam exponential family.
Fungsi Link
Fungsi link kanonik untuk distribusi Poisson adalah fungsi logaritma:
\[ g(\mu) = \log(\mu) \]
Sehingga modelnya menjadi:
\[ \log(\mu_i) = x_i^T \beta \]
dan fungsi inverse link:
\[ \mu_i = \exp(x_i^T \beta) \]
Estimasi Parameter
Estimasi parameter \(( \beta )\)
dilakukan dengan metode Maximum Likelihood Estimation
(MLE).
Log-likelihood fungsi untuk regresi Poisson:
\[ l(\beta) = \sum_{i=1}^n \left[ y_i x_i^T \beta - \exp(x_i^T \beta) - \log(y_i!) \right] \]
Nilai \(( \beta )\) dapat diperoleh melalui metode numerik seperti iterasi Newton-Raphson.
Contoh
Misalkan kita memiliki data sebagai berikut:
# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)
# Buat data biner
game <- rep(1:23, Attempts)
result <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))
x <- data$Game
y <- data$Result
# Gabung jadi data frame
data <- data.frame(x, y)
head(data)
lambda <- exp(0.3 + 0.6 * x)
Estimasi Regresi Poisson
poisson_model <- glm(y ~ x, data = data, family = poisson())
summary(poisson_model)
##
## Call:
## glm(formula = y ~ x, family = poisson(), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.679572 0.189073 -3.594 0.000325 ***
## x -0.008375 0.013544 -0.618 0.536344
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 211.97 on 295 degrees of freedom
## Residual deviance: 211.59 on 294 degrees of freedom
## AIC: 485.59
##
## Number of Fisher Scoring iterations: 5
Plot Hasil Prediksi
plot(x, y, pch = 16, col = "darkgray", main = "Data dan Hasil Prediksi")
newdata <- data.frame(x=seq(min(x), max(x), length.out = 100))
pred <- predict(poisson_model, newdata = newdata, type = "response")
lines(newdata$x, pred, col = "blue", lwd = 2)
Diagnostik dan Overdispersion
Salah satu asumsi penting dari model Poisson adalah bahwa mean dan varians dari variabel respons adalah sama (𝐸[𝑌]=𝑉𝑎𝑟[𝑌]). Jika varians lebih besar dari mean, maka terjadi overdispersion. Untuk mendeteksi overdispersion:
dispersion <- sum(residuals(poisson_model, type = "pearson")^2) / poisson_model$df.residual
dispersion
## [1] 0.54753
Jika nilai dispersion > 1, maka overdispersion mungkin terjadi dan model alternatif seperti Negative Binomial Regression dapat digunakan.
exp(coef(poisson_model))
## (Intercept) x
## 0.5068336 0.9916599
Visualisasi Prediksi
data$predicted <- predict(poisson_model, type = "response")
library(ggplot2)
ggplot(data, aes(x=x, y=y)) +
geom_jitter(width = 0.2, alpha = 0.6) +
geom_point(aes(y=predicted), shape = 18, size = 3, color = "black") +
labs(title = "Prediksi", x="Game", y="Result") +
theme_minimal()
Evaluasi Model
# Plot residuals
plot(poisson_model$residuals,
main = "Residual Plot",
ylab = "Residual",
xlab = "Index",
pch = 19,
col = "steelblue")
abline(h = 0, col = "red", lty = 2)
Dalam Generalized Linear Model (GLM), inferensi statistik membutuhkan pemahaman terhadap ekspektasi dan variansi dari estimator model, terutama untuk mengembangkan alat-alat uji seperti Wald test, Likelihood Ratio test, dan interval kepercayaan.
Ekspektasi dan Varians dalam GLM
Ekspektasi menunjukkan apakah suatu estimator tak bias, yaitu:
\[ \mathbb{E}[\hat{\beta}] = \beta \]
Dalam GLM, MLE dari \(\hat{\beta}\) bersifat asymptotically unbiased.
Varians menunjukkan presisi dari estimasi parameter:
\[ \text{Var}(\hat{\beta}) \approx (\mathbf{X}^T \mathbf{W} \mathbf{X})^{-1} \]
dimana \(\mathbf{W}\) adalah matriks bobot yang tergantung pada distribusi dan fungsi link.
Distribusi Asimptotik Estimator
Dengan ukuran sampel besar:
\[ \hat{\beta} \sim \mathcal{N}(\beta, \text{Var}(\hat{\beta})) \]
Distribusi ini adalah dasar dari:
Varians dalam GLM Tidak Konstan
Tidak seperti regresi linear (OLS) yang mengasumsikan homoskedastisitas:
\[ \text{Var}(Y_i) = \sigma^2 \]
Dalam GLM:
\[ \text{Var}(Y_i) = \phi V(\mu_i) \]
Contoh:
Contoh Regresi Poisson
# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)
# Buat data biner
x <- rep(1:23, Attempts)
y <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))
# Gabung jadi data frame
data <- data.frame(x = x, y = y)
# Ringkasan model regresi Poisson
model <- glm(y ~ x, data = data, family = poisson())
summary(model)
##
## Call:
## glm(formula = y ~ x, family = poisson(), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.679572 0.189073 -3.594 0.000325 ***
## x -0.008375 0.013544 -0.618 0.536344
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 211.97 on 295 degrees of freedom
## Residual deviance: 211.59 on 294 degrees of freedom
## AIC: 485.59
##
## Number of Fisher Scoring iterations: 5
Kesimpulan:
Ekspektasi
Jika diturunkan berdasarkan fungsi momen:
\[ E(Y) = \int y f(y; \theta) \, dy = \mu \]
Untuk keluarga eksponensial:
\[ \log f(y; \theta) = a(y) + b(\theta) y + c(\theta) \]
atau:
\[ \log f(y; \theta) = y \theta - b(\theta) + c(y) \]
Maka ekspektasi turunan pertama:
\[ U(\theta) = \frac{\partial \ell}{\partial \theta} = y - b'(\theta) \] Dan ekspektasi turunan pertama:
\[ \mathbb{E}[U(\theta)] = \mathbb{E}[y - b'(\theta)] = \mu - b'(\theta) = 0 \]
Maka:
\[ \mu = b'(\theta) \]
Varians
Turunan kedua:
\[ \frac{\partial^2 \ell}{\partial \theta^2} = -b''(\theta) \]
Sehingga:
\[ \text{Var}(Y) = b''(\theta) = \phi V(\mu) \]
Maximum Likelihood Estimation (MLE)
Namun, karena bentuk GLM tidak eksplisit, digunakan metode numerik.
Metode Optimisasi: Newton-Raphson
Iterasi:
\[ \beta^{(t+1)} = \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \] Fisher Scoring
IRLS (Iteratively Reweighted Least Square) - Modifikasi dari Fisher scoring, hasil estimasi mirip dengan Least Square.
Implementasi Newton-Raphson Statistik score ke-𝑗 Statistik score ke-𝑗:
\[ U_j(\beta) = \frac{\partial \log L(\beta)}{\partial \beta_j} \] Turunan kedua:
\[ H_{jk}(\beta) = \frac{\partial^2 \log L(\beta)}{\partial \beta_j \partial \beta_k} \]
Taylor expansion:
\[ U(\beta^*) \approx U(\beta) + H(\beta)(\beta^* - \beta) \]
Estimasi parameter:
\[ \hat{\beta} \approx \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \]
Diagnostik digunakan untuk mengevaluasi apakah model sudah tepat.
Statistik Devians
\[ D = 2 \sum \left[ y_i \log \left( \frac{y_i}{\hat{\mu}_i} \right) - (y_i - \hat{\mu}_i) \right] \]
Statistik Chi-Kuadrat Pearson
\[ X^2 = \sum \left( \frac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i} \right) \] • Jika signifikan → model lebih baik daripada tanpa model.
Catatan
• Untuk data yang dikelompokkan, statistik devians dan chi-kuadrat Pearson mengikuti distribusi ChiSquare.
• Untuk data tidak dikelompokkan, tidak mengikuti distribusi Chi-Square.
• Devians diminimalkan oleh MLE → cocok digunakan untuk evaluasi model.
Analisis Residual
• Residual adalah selisih antara observasi dengan prediksi.
• Dapat digunakan untuk memeriksa penyimpangan sistematis.
• Dapat diplot untuk menilai asumsi model
Regresi logistik digunakan untuk memodelkan probabilitas dari variabel respons biner (0/1) berdasarkan satu atau lebih variabel prediktor. Estimasi parameter dilakukan menggunakan Maximum Likelihood Estimation (MLE) karena model tidak linear dalam parameternya.
Fungsi model logistik:
\[ \pi(x) = \frac{\exp(\beta_0 + \beta_1 x)}{1 + \exp(\beta_0 + \beta_1 x)} \]
Log-likelihood untuk \(( n )\) observasi:
\[ \ell(\beta) = \sum_{i=1}^{n} \left[ y_i \log(\pi_i) + (1 - y_i) \log(1 - \pi_i) \right] \] Estimasi dengan Newton-Raphson
Metode Newton-Raphson digunakan untuk mencari nilai parameter \(( \beta )\) yang memaksimalkan fungsi log-likelihood pada model regresi logistik.
Fungsi Log-Likelihood
Model regresi logistik untuk probabilitas:
\[ \pi_i = \frac{1}{1 + \exp(-x_i^\top \beta)} \]
Log-likelihood untuk \(( n )\) observasi:
\[ \ell(\beta) = \sum_{i=1}^{n} \left[ y_i \log(\pi_i) + (1 - y_i) \log(1 - \pi_i) \right] \]
Langkah-langkah Newton-Raphson Turunan dalam Estimasi Logistik
\[ U(\beta) = \frac{\partial \ell(\beta)}{\partial \beta} = X^\top (y - \pi) \]
\[ H(\beta) = -X^\top W X, \quad \text{dengan } W = \text{diag}(\pi_i (1 - \pi_i)) \]
\[ \beta^{(t+1)} = \beta^{(t)} + (X^\top W^{(t)} X)^{-1} X^\top (y - \pi^{(t)}) \]
Estimasi MLE dengan Newton-Raphson (Manual di R)
# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)
# Buat data biner
x <- rep(1:23, Attempts)
y <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))
# Gabung jadi data frame
data <- data.frame(x = x, y = y)
beta_true <- c(-1, 2)
X <- cbind(1, x)
eta <- X %*% beta_true
p <- 1 / (1 + exp(-eta))
Newton-Raphson Iterasi Manual
beta <- c(0, 0)
tol <- 1e-6
max_iter <- 100
for (i in 1:max_iter) {
eta <- X %*% beta
pi_hat <- 1 / (1 + exp(-eta))
W <- diag(as.numeric(pi_hat * (1 - pi_hat)))
z <- eta + solve(W) %*% (y - pi_hat)
beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z
if (sum(abs(beta_new - beta)) < tol) break
beta <- beta_new }
beta # hasil estimasi akhir
## [,1]
## 0.02085078
## x -0.01546112
• Estimasi parameter pada model regresi logistik dilakukan dengan MLE. • Newton-Raphson adalah metode numerik yang digunakan untuk memaksimalkan log-likelihood.
• Iterasi didasarkan pada turunan pertama (score) dan kedua (Hessian).
• Prosedur identik dengan IRLS (Iteratively Reweighted Least Squares) dalam implementasi GLM.
Inferensi Parameter
Tujuan Uji Wald
Untuk menguji signifikansi parameter \((
\beta_j )\) dalam model regresi logistik:
Teori Wald Test
Dari teori estimasi MLE, estimator \((
\hat{\beta}_j )\) mendekati distribusi normal:
\[ \hat{\beta}_j \sim N(\beta_j, \text{Var}(\hat{\beta}_j)) \]
Jika \(( H_0 )\) benar (yaitu \(( \beta_j = 0 )\)), maka:
\[ Z = \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)} \sim N(0,1) \]
Dengan Statistik Wald:
\[ W = Z^2 = \left( \frac{\hat{\beta}_j}{SE(\hat{\beta}_j)} \right)^2 \sim \chi^2_1 \]
Uji Wald Langkah demi Langkah
# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)
# Buat data biner
x <- rep(1:23, Attempts)
y <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))
# Gabung jadi data frame
data <- data.frame(x = x, y = y)
log_odds <- -0.5 + 1.2 * x
p <- 1 / (1 + exp(-log_odds))
model <- glm(y ~ x, data = data, family = binomial)
summary(model)
##
## Call:
## glm(formula = y ~ x, family = binomial, data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.02085 0.26206 0.080 0.937
## x -0.01546 0.01845 -0.838 0.402
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 408.06 on 295 degrees of freedom
## Residual deviance: 407.35 on 294 degrees of freedom
## AIC: 411.35
##
## Number of Fisher Scoring iterations: 3
Langkah 1: Ambil nilai koefisien dan standard error (SE)
beta_hat <- coef(model)["x"]
se_beta <- summary(model)$coefficients["x", "Std. Error"]
Langkah 2: Hitung Statistik Z
Z <- beta_hat / se_beta
Z
## x
## -0.8380593
Langkah 3: Hitung Statistik Wald
Wald_stat <- Z^2
Wald_stat
## x
## 0.7023434
Langkah 4: Hitung p-value
p_value <- 1 - pchisq(Wald_stat, df = 1)
p_value
## x
## 0.4019974
Interpretasi
• Jika p-value < 0.05, maka koefisien signifikan → variabel prediktor berpengaruh.
• Jika p-value > 0.05, maka tidak ada cukup bukti untuk menolak - \(( H_0 )\)
Kesimpulan
Uji Wald didasarkan pada rasio antara estimasi parameter dan standar error-nya. Dengan menaikkan nilai Z menjadi kuadrat (Z²), kita memperoleh distribusi Chi-Square untuk pengujian hipotesis parameter individual dalam model regresi logistik
Bandingkan model penuh dengan model tanpa prediktor
# Model null
model_null <- glm(y ~ 1, data = data, family = binomial)
# Likelihood ratio test
anova(model_null, model, test = "Chisq")
Evaluasi Kebaikan Model
Semakin kecil AIC, semakin baik model.
AIC(model)
## [1] 411.3528
Alternatif terhadap AIC, menghukum kompleksitas model.
BIC(model)
## [1] 418.7335
• Estimasi regresi logistik dilakukan dengan MLE melalui iterasi Newton-Raphson.
• Inferensi parameter dapat dilakukan dengan uji Wald dan likelihood ratio (uji Chi-Square).
• AIC dan BIC digunakan untuk mengevaluasi kompleksitas dan kecocokan model.
Model regresi Poisson digunakan untuk memodelkan data count (jumlah kejadian) dimana variabel respons mengikuti distribusi Poisson. Estimasi dilakukan dengan Maximum Likelihood Estimation (MLE), dan inferensi dilakukan dengan uji Wald dan Likelihood Ratio Test.
Model Regresi Poisson
Distribusi Poisson:
\[ P(Y_i = y_i) = \frac{e^{-\lambda_i} \lambda_i^{y_i}}{y_i!} \]
Model regresi Poisson:
\[ \log(\lambda_i) = x_i^T \beta \] Estimasi Parameter (MLE)
Log-likelihood fungsi:
\[ \ell(\beta) = \sum_{i=1}^{n} \left[ y_i \log(\lambda_i) - \lambda_i - \log(y_i!) \right] \]
Dengan:
\[ \lambda_i = \exp(x_i^T \beta) \]
Estimasi dilakukan dengan metode iterasi (IRLS)
Estimasi parameter model regresi Poisson menggunakan pendekatan
Maximum Likelihood Estimation (MLE) dengan metode
Iteratively Reweighted Least Squares (IRLS) secara
manual, tanpa menggunakan glm()
.
Tahap 1: Definisikan Model Regresi Poisson
\[ \log(\lambda_i) = x_i^T \beta \quad \Rightarrow \quad \lambda_i = \exp(x_i^T \beta) \]
Tahap 2: Mencari Log-likelihood yang dimaksimumkan
\[ \ell(\beta) = \sum_{i=1}^{n} \left[ y_i \log(\lambda_i) - \lambda_i - \log(y_i!) \right] \] Tahap 3: Formulasi iteratif:
\[ \beta^{(t+1)} = \left( \mathbf{X}^T \mathbf{W}^{(t)} \mathbf{X} \right)^{-1} \mathbf{X}^T \mathbf{W}^{(t)} \mathbf{z}^{(t)} \]
Dengan:
dan
\[ \eta_i = \log(\lambda_i) = x_i^T \beta \] Simulasi Data
# Data agregat
Made <- c(4, 5, 5, 5, 2, 7, 6, 9, 4, 1, 13, 5, 6, 9, 7, 3, 8, 1, 18, 3, 10, 1, 3)
Attempts <- c(5, 11, 14, 12, 7, 10, 14, 15, 12, 4, 27, 17, 12, 9, 12, 10, 12, 6, 39, 13, 17, 6, 12)
# Buat data biner
x <- rep(1:23, Attempts)
y <- unlist(mapply(function(m, a) c(rep(1, m), rep(0, a - m)), Made, Attempts))
# Gabung jadi data frame
data <- data.frame(x = x, y = y)
X <- cbind(1, x) # Tambah intercept
beta_true <- c(0.5, 0.8)
eta <- X %*% beta_true
lambda <- exp(eta)
IRLS Manual Step-by-Step
# Inisialisasi
beta <- c(0, 0)
tol <- 1e-6
max_iter <- 100
for (i in 1:max_iter) {
eta <- X %*% beta
lambda <- exp(eta)
W <- diag(as.numeric(lambda))
z <- eta + (y - lambda) / lambda
beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z
if (sum(abs(beta_new - beta)) < tol) {
cat("Konvergen pada iterasi ke-", i, "\n")
break
}
beta <- beta_new
}
## Konvergen pada iterasi ke- 5
beta # hasil estimasi
## [,1]
## -0.679572484
## x -0.008375021
Perbandingan dengan glm
model_glm <- glm(y ~ x, family = poisson())
summary(model_glm)
##
## Call:
## glm(formula = y ~ x, family = poisson())
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.679572 0.189073 -3.594 0.000325 ***
## x -0.008375 0.013544 -0.618 0.536344
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 211.97 on 295 degrees of freedom
## Residual deviance: 211.59 on 294 degrees of freedom
## AIC: 485.59
##
## Number of Fisher Scoring iterations: 5
• IRLS memberikan cara iteratif untuk menghitung estimasi MLE dalam regresi Poisson.
• Hasil manual IRLS sangat mendekati hasil glm() dari R.
• Metode ini memberikan pemahaman mendalam atas mekanisme di balik fungsi glm().
Pengujian hipotesis Uji Wald
Untuk menguji H0: = 0
# Koefisien dan standar error
coef_val <- coef(model)[2]
se_val <- summary(model)$coefficients[2, 2]
wald_z <- coef_val / se_val
wald_chisq <- wald_z^2
p_value <- 1 - pchisq(wald_chisq, df = 1)
cat("Z:", wald_z, "\nChi-Square:", wald_chisq, "\np-value:", p_value)
## Z: -0.8380593
## Chi-Square: 0.7023434
## p-value: 0.4019974
Uji Likelihood Ratio (Chi-Square)
model_null <- glm(y ~ 1, family = poisson(), data = data)
anova(model_null, model, test = "Chisq")
Evaluasi Model (AIC & BIC)
AIC(model)
## [1] 411.3528
BIC(model)
## [1] 418.7335
• Estimasi parameter regresi Poisson dilakukan menggunakan MLE
• Uji Wald dan Likelihood Ratio digunakan untuk pengujian hipotesis
• AIC dan BIC digunakan untuk evaluasi dan pemilihan model terbaik
Regresi logistik merupakan salah satu metode analisis statistik yang digunakan untuk memodelkan hubungan antara satu variabel respons biner (dua kategori) dengan satu atau lebih variabel prediktor. Dalam penerapannya, prediktor yang digunakan dapat memiliki skala pengukuran berbeda, yaitu:
Nominal: Variabel yang tidak memiliki urutan logis antar kategorinya, hanya sebagai pembeda.
Contoh: jenis kelamin (laki-laki, perempuan). Dalam analisis regresi logistik, data nominal diubah menjadi variabel dummy (misal: Female = 1, Male = 0).
Ordinal: Variabel yang memiliki urutan logis antar kategori, tetapi jarak antar kategori belum tentu sama.
Contoh: tingkat pendidikan (SMA < Sarjana < Master < Doktor). Dalam analisis, variabel ordinal bisa diperlakukan:
– Sebagai nominal dengan membuat dummy untuk setiap kategori.
– Seperti rasio dengan memberi nilai tertentu (misalnya SMA = 1, Sarjana = 2, Master = 3, Doktor = 4) untuk mencerminkan urutan dan jarak antar tingkatannya.
Rasio: Variabel numerik kontinu yang memiliki nol absolut dan rasio bermakna.
Contoh: pendapatan (dalam juta rupiah per bulan). Data ini dapat langsung digunakan dalam model tanpa transformasi khusus.
Dalam dunia bisnis modern, startup merupakan salah satu bentuk usaha yang paling dinamis namun juga paling rentan mengalami kegagalan. Untuk memahami faktor-faktor yang memengaruhi keberhasilan atau kegagalan sebuah startup, dilakukan simulasi data terhadap 250 proyek startup dari berbagai bidang. Setiap proyek dikategorikan berdasarkan bidang usaha (field), tingkat pengalaman founder (experience), serta jumlah modal awal (capital). Tujuan analisis ini adalah untuk memodelkan probabilitas kegagalan proyek menggunakan regresi logistik biner, di mana kegagalan didefinisikan sebagai 1 dan keberhasilan sebagai 0.
Variabel-variabel yang digunakan dalam studi ini adalah:
Failure (biner): Status hasil proyek startup. Nilai 1 menunjukkan startup gagal, sedangkan nilai 0 menunjukkan startup berhasil.
Field (kategori): Bidang usaha startup, terdiri atas empat kategori, yaitu Fintech, Healthtech, Edtech, dan E-commerce.
Experience (kategori): Tingkat pengalaman pendiri (founder) startup, dikategorikan menjadi Low, Medium, dan High.
Capital (numerik): Modal awal startup yang dimiliki pada saat pendirian, dinyatakan dalam ribuan dolar Amerika Serikat (USD).
Dengan menggunakan ketiga variabel prediktor tersebut, probabilitas kegagalan setiap proyek dianalisis untuk melihat bagaimana bidang usaha, pengalaman pendiri, dan besar modal awal berkontribusi terhadap risiko kegagalan startup
Simulasi dataset akan dibuat menggunakan 250 observasi
library(tibble)
library(dplyr)
set.seed(123)
# Simulasi variabel
n <- 250
# Variabel X: bidang startup
field <- sample(c("Fintech", "Healthtech", "Edtech", "E-commerce"),
n, replace = TRUE, prob = c(0.3, 0.2, 0.2, 0.3))
# Variabel Y: tingkat pengalaman founder
experience <- sample(c("Low", "Medium", "High"),
n, replace = TRUE, prob = c(0.4, 0.4, 0.2))
# Variabel tambahan: jumlah modal awal (dalam ribu USD)
capital <- rnorm(n, mean = 500, sd = 100)
# Membuat peluang kegagalan
logit_p <- -1.5 +
0.7 * (field == "E-commerce") +
0.5 * (experience == "Low") +
0.002 * (500 - capital) # semakin kecil modal, semakin besar peluang gagal
p <- 1 / (1 + exp(-logit_p))
set.seed(123)
failure <- rbinom(n, 1, p) # 1 = gagal, 0 = berhasil
# Gabungkan menjadi data frame
library(tibble)
startup_data <- tibble(failure, field, experience, capital)
head(startup_data)
startup_data %>%
dplyr::group_by(failure) %>%
dplyr::summarise(
Jumlah = dplyr::n(),
Rata2_Modal = mean(capital)
)
startup_data_nominal <- startup_data %>%
mutate(
field = factor(field, levels = c("Fintech", "Healthtech", "Edtech", "E-commerce")),
experience = factor(experience, levels = c("Low", "Medium", "High"))
)
model_nominal <- glm(failure ~ field + experience + capital,
data = startup_data_nominal,
family = binomial)
summary(model_nominal)
##
## Call:
## glm(formula = failure ~ field + experience + capital, family = binomial,
## data = startup_data_nominal)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.260e+01 1.848e+03 -0.007 0.99456
## fieldHealthtech 2.659e+01 1.848e+03 0.014 0.98852
## fieldEdtech 1.824e+01 1.848e+03 0.010 0.99212
## fieldE-commerce 1.758e+01 1.848e+03 0.010 0.99241
## experienceMedium -3.280e+00 1.157e+00 -2.835 0.00458 **
## experienceHigh -5.289e+00 1.703e+00 -3.107 0.00189 **
## capital -1.422e-02 4.432e-03 -3.208 0.00134 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 260.899 on 249 degrees of freedom
## Residual deviance: 59.771 on 243 degrees of freedom
## AIC: 73.771
##
## Number of Fisher Scoring iterations: 19
Interpretasi Koefisien
Intercept (-12.60)
-Ini adalah log-odds dasar untuk startup yang bergerak di bidang Fintech, dengan founder berpengalaman Low, dan modal standar (capital = 0).
-Nilai p = 0.99456 (tidak signifikan di level 5%), artinya baseline ini tidak berbeda signifikan terhadap probabilitas gagal.
-Odds Ratio = exp(-12.60) ≈ 0.00000336 → Sangat kecil. Peluang kegagalan pada kondisi baseline hampir nol, meski tidak signifikan.
fieldHealthtech (+26.59)
-Startup di bidang Healthtech memiliki log-odds kegagalan lebih tinggi sebesar 26.59 dibandingkan Fintech (referensi).
-p = 0.98852 (tidak signifikan), artinya tidak ada bukti signifikan bahwa Healthtech lebih berisiko gagal dibanding Fintech.
-Odds Ratio = exp(26.59) ≈ 3.46 × 10¹¹ → Secara matematis sangat besar, tapi tidak signifikan.
fieldEdtech (+18.24)
-Startup di bidang Edtech memiliki log-odds kegagalan lebih tinggi sebesar 18.24 dibandingkan Fintech.
-p = 0.99212 (tidak signifikan), sehingga perbedaan tidak cukup bukti untuk menyatakan risiko kegagalan Edtech lebih tinggi dibanding Fintech.
-Odds Ratio = exp(18.24) ≈ 9.2 × 10⁷ → Besar, tapi tidak signifikan.
fieldE-commerce (+17.58)
-Startup di bidang E-commerce memiliki log-odds kegagalan lebih tinggi sebesar 17.58 dibandingkan Fintech.
-p = 0.99241 (tidak signifikan), sehingga tidak cukup bukti bahwa E-commerce lebih rentan gagal dibanding Fintech.
-Odds Ratio = exp(17.58) ≈ 4.3 × 10⁷ → Besar, namun tidak signifikan.
experienceMedium (-3.28)
-Founder dengan pengalaman Medium memiliki log-odds kegagalan lebih rendah sebesar 3.28 dibandingkan Low (referensi).
-p = 0.00458 (signifikan di level 5%), berarti pengalaman Medium signifikan menurunkan peluang kegagalan.
-Odds Ratio = exp(-3.28) ≈ 0.038 → Peluang gagal hanya sekitar 3.8% dibandingkan dengan founder berpengalaman Low.
experienceHigh (-5.289)
-Founder dengan pengalaman High memiliki log-odds kegagalan lebih rendah sebesar 5.289 dibandingkan Low.
-p = 0.00189 (signifikan di level 5%), menunjukkan bahwa pengalaman tinggi sangat signifikan mengurangi kegagalan.
-Odds Ratio = exp(-5.289) ≈ 0.005 → Peluang gagal hanya sekitar 0.5% dari peluang founder berpengalaman rendah.
capital (-0.01422)
-Setiap kenaikan modal sebesar 1 unit (asumsinya dalam satuan modal standar), menurunkan log-odds kegagalan sebesar 0.01422.
-p = 0.00134 (signifikan di level 5%), artinya modal lebih besar signifikan menurunkan peluang kegagalan.
-Odds Ratio = exp(-0.01422) ≈ 0.9859 → Setiap tambahan modal 1 unit mengurangi peluang gagal sekitar 1.4%.
Interpretasi Goodness-of-Fit
Null Deviance (260.899): Deviance model tanpa prediktor.
Residual Deviance (59.771): Deviance model dengan prediktor. → Penurunan dari null deviance ke residual deviance menunjukkan bahwa model membawa informasi yang cukup baik, meski penurunannya sedang saja (tidak drastis).
AIC (73.771): → Digunakan untuk membandingkan model: semakin kecil AIC, semakin baik model dalam menyeimbangkan goodness-of-fit dan kompleksitas. → Nilai AIC ini tergolong cukup rendah untuk data simulasi seperti ini.
Signifikansi Model
Variabel yang signifikan: → experienceMedium, experienceHigh, dan capital berpengaruh signifikan terhadap penurunan peluang gagal startup (p < 0.01).
Variabel yang tidak signifikan: → Variabel field (baik Healthtech, Edtech, maupun E-commerce) tidak signifikan mempengaruhi peluang kegagalan dibandingkan Fintech.
Arah hubungan: → Semakin tinggi pengalaman dan modal, semakin kecil peluang startup untuk mengalami kegagalan.
Kesimpulan Praktis
Pengalaman founder (medium atau high) merupakan faktor kunci dalam menurunkan risiko kegagalan startup.
Modal yang lebih besar memberikan perlindungan signifikan terhadap kemungkinan gagal.
Bidang usaha (field) startup tidak menunjukkan pengaruh yang signifikan terhadap kegagalan dalam model ini, menunjukkan bahwa faktor internal (seperti pengalaman dan modal) lebih berperan.
Model sudah lebih baik dibandingkan model null, tapi masih ada ruang untuk peningkatan model, misalnya dengan menambahkan variabel baru atau memperbaiki struktur data.
startup_data_numeric <- startup_data %>%
mutate(
experience_numeric = case_when(
experience == "Low" ~ 1,
experience == "Medium" ~ 2,
experience == "High" ~ 3
)
)
model_numeric <- glm(failure ~ field + experience_numeric + capital,
data = startup_data_numeric,
family = binomial)
summary(model_numeric)
##
## Call:
## glm(formula = failure ~ field + experience_numeric + capital,
## family = binomial, data = startup_data_numeric)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 7.684e+00 2.372e+00 3.240 0.001197 **
## fieldEdtech 6.730e-01 8.064e-01 0.835 0.403978
## fieldFintech -1.758e+01 1.861e+03 -0.009 0.992466
## fieldHealthtech 9.160e+00 1.753e+00 5.224 1.75e-07 ***
## experience_numeric -2.799e+00 8.325e-01 -3.362 0.000773 ***
## capital -1.417e-02 4.363e-03 -3.247 0.001168 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 260.899 on 249 degrees of freedom
## Residual deviance: 60.258 on 244 degrees of freedom
## AIC: 72.258
##
## Number of Fisher Scoring iterations: 19
Interpretasi Koefisien
Intercept (+7.684)
Ini adalah log-odds dasar untuk startup pada bidang E-commerce (sebagai baseline), dengan pengalaman dan modal nol. Nilai ini signifikan (p = 0.0012), menunjukkan bahwa probabilitas kegagalan dasar cukup tinggi.
fieldEdtech (+0.6730)
Startup pada bidang Edtech memiliki log-odds kegagalan yang lebih tinggi sebesar 0.673 dibandingkan E-commerce, namun tidak signifikan (p = 0.4040). Ini menunjukkan bahwa perbedaan antara Edtech dan E-commerce dalam hal kegagalan tidak cukup kuat untuk disimpulkan secara statistik.
fieldFintech (–17.5800)
Startup pada bidang Fintech memiliki log-odds kegagalan yang jauh lebih rendah dibandingkan E-commerce, namun tidak signifikan (p = 0.9925). Nilai koefisien ekstrem ini kemungkinan muncul akibat data yang sangat tidak seimbang (jumlah Fintech sangat sedikit atau tidak gagal).
fieldHealthtech (+9.1600)
Startup pada bidang Healthtech memiliki log-odds kegagalan lebih tinggi secara signifikan sebesar 9.16 dibandingkan E-commerce. p < 0.001, sangat signifikan.
Odds Ratio ≈ exp(9.16) ≈ 9527.55 Artinya, peluang gagal startup Healthtech sekitar 9527 kali lebih tinggi dari E-commerce (namun interpretasi ini harus hati-hati karena nilainya sangat ekstrem dan bisa dipengaruhi data outlier/struktur data tertentu).
experience_numeric (–2.799)
Setiap kenaikan 1 tingkat pengalaman (misalnya dari Low → Medium, atau Medium → High) menurunkan log-odds kegagalan sebesar 2.799.
p = 0.0008 (sangat signifikan)
Odds Ratio = exp(–2.799) ≈ 0.061 Artinya, setiap peningkatan pengalaman menurunkan peluang gagal sekitar 94%.
capital (–0.01417)
Setiap penambahan 1 unit modal (misalnya 1 juta), menurunkan log-odds kegagalan sebesar 0.0142.
p = 0.0012 (sangat signifikan)
Odds Ratio = exp(–0.01417) ≈ 0.986 Artinya, setiap tambahan 1 juta modal menurunkan peluang gagal sekitar 1.4%.
Interpretasi Goodness-of-Fit
Null deviance (260.899): Deviance model tanpa prediktor.
Residual deviance (60.258): Deviance model dengan prediktor. → Penurunan dari null deviance ke residual deviance menunjukkan bahwa prediktor memberikan tambahan informasi yang sangat berarti.
AIC (72.258) → Digunakan untuk membandingkan model: semakin kecil AIC, semakin baik keseimbangan antara goodness-of-fit dan kompleksitas model.
Signifikansi Model
Variabel experience_numeric (p = 0.000773) dan capital (p = 0.001168) sangat signifikan dalam menurunkan peluang kegagalan.
Variabel fieldHealthtech (p = 1.75e-07) juga signifikan, dengan log-odds kegagalan jauh lebih tinggi dibandingkan bidang referensi (E-commerce).
Variabel fieldEdtech dan fieldFintech tidak signifikan secara statistik, artinya tidak terbukti berbeda secara nyata dari bidang E-commerce.
Kesimpulan Praktis
Semakin tinggi tingkat pengalaman, semakin kecil peluang startup mengalami kegagalan.
Modal (capital) yang lebih besar secara konsisten menurunkan peluang kegagalan, walaupun efeknya kecil (koefisien -0.014).
Startup di bidang Healthtech memiliki risiko kegagalan jauh lebih tinggi dibandingkan bidang E-commerce (log-odds = 9.16).
Variabel bidang lain seperti Edtech dan Fintech tidak menunjukkan perbedaan signifikan dalam hal peluang kegagalan.
Model ini menunjukkan performa yang sangat baik dalam menjelaskan kegagalan startup, dan masih bisa ditingkatkan dengan menambah variabel lain jika tersedia
Perbandingan Model
list(
AIC_Nominal = AIC(model_nominal),
AIC_Numeric = AIC(model_numeric)
)
## $AIC_Nominal
## [1] 73.77055
##
## $AIC_Numeric
## [1] 72.25842
Keterangan: Semakin kecil nilai AIC, semakin bagus untuk dipakai
Goodness-of-Fit Model
nullmod <- glm(failure ~ 1, data = startup_data, family = binomial)
r2_nominal <- 1 - (logLik(model_nominal)/logLik(nullmod))
r2_numeric <- 1 - (logLik(model_numeric)/logLik(nullmod))
list(
McFadden_R2_Nominal = r2_nominal,
McFadden_R2_Numeric = r2_numeric
)
## $McFadden_R2_Nominal
## 'log Lik.' 0.7709056 (df=7)
##
## $McFadden_R2_Numeric
## 'log Lik.' 0.7690357 (df=6)
Interpretasi: Semakin besar nilainya, semakin baik model memprediksi dibandingkan model null
Visualisasi Prediksi
library(ggplot2)
startup_data_nominal <- startup_data_nominal %>%
mutate(predicted = predict(model_nominal, type = "response"))
startup_data_numeric <- startup_data_numeric %>%
mutate(predicted = predict(model_numeric, type = "response"))
# Plot untuk model nominal
startup_data_nominal %>%
ggplot(aes(x = capital, y = predicted, color = experience)) +
geom_point(alpha = 0.6) +
labs(
title = "Prediksi Probabilitas (Ordinal sebagai Nominal)",
x = "Capital",
y = "Prediksi Probabilitas"
) +
theme_minimal()
# Plot untuk model numeric
startup_data_numeric %>%
ggplot(aes(x = capital, y = predicted, color = experience)) +
geom_point(alpha = 0.6) +
labs(
title = "Prediksi Probabilitas (Ordinal sebagai Numeric)",
x = "Capital",
y = "Prediksi Probabilitas"
) +
theme_minimal()
Ringkasan Koefisien Model
term | estimate | std.error | statistic | p.value | conf.low | conf.high |
---|---|---|---|---|---|---|
(Intercept) | -12.600 | 1847.902 | -0.007 | 0.995 | NA | 118.728 |
fieldHealthtech | 26.591 | 1847.901 | 0.014 | 0.989 | -38.652 | 710.624 |
fieldEdtech | 18.240 | 1847.901 | 0.010 | 0.992 | -113.088 | NA |
fieldE-commerce | 17.581 | 1847.901 | 0.010 | 0.992 | -113.747 | NA |
experienceMedium | -3.280 | 1.157 | -2.835 | 0.005 | -6.323 | -1.395 |
experienceHigh | -5.289 | 1.703 | -3.107 | 0.002 | -9.349 | -2.472 |
capital | -0.014 | 0.004 | -3.208 | 0.001 | -0.024 | -0.006 |
term | estimate | std.error | statistic | p.value | conf.low | conf.high |
---|---|---|---|---|---|---|
(Intercept) | 7.684 | 2.372 | 3.240 | 0.001 | 3.441 | 12.907 |
fieldEdtech | 0.673 | 0.806 | 0.835 | 0.404 | -0.949 | 2.290 |
fieldFintech | -17.575 | 1861.204 | -0.009 | 0.992 | NA | 114.728 |
fieldHealthtech | 9.160 | 1.753 | 5.224 | 0.000 | 6.323 | 13.388 |
experience_numeric | -2.799 | 0.833 | -3.362 | 0.001 | -4.803 | -1.425 |
capital | -0.014 | 0.004 | -3.247 | 0.001 | -0.024 | -0.006 |
Kesimpulan
-Semakin tinggi tingkat pengalaman, semakin kecil peluang startup mengalami kegagalan.
-Modal yang lebih besar memberikan perlindungan signifikan terhadap kemungkinan gagal.
-Bidang usaha (field) startup tidak menunjukkan pengaruh yang signifikan terhadap kegagalan dalam model ini, menunjukkan bahwa faktor internal (seperti pengalaman dan modal) lebih berperan.
Dalam analisis regresi logistik, pemilihan model sangat krusial untuk mendapatkan model yang baik dalam memprediksi probabilitas kejadian suatu peristiwa (respon biner). Dua pendekatan utama dalam membangun model adalah pendekatan Confirmatory dan Exploratory
1. Confirmatory (Pendekatan Konfirmator)
Pendekatan ini digunakan ketika peneliti telah memiliki teori atau hipotesis yang jelas mengenai efek atau hubungan antara variabel prediktor dan respon.
Ciri-ciri:
• Model dibangun berdasarkan teori atau hasil penelitian sebelumnya.
• Tujuan utamanya adalah menguji apakah efek tersebut benar-benar signifkan, bukan sekadar mencari model terbaik.
• Peneliti biasanya menyusun model penuh terlebih dahulu, lalu menguji apakah penambahan atau pengurangan suatu efek (misalnya, interaksi) memberikan peningkatan model secara signifkan.
• Uji signifkansi dilakukan dengan membandingkan model dengan efek tertentu dan model tanpa efek tersebut, misalnya dengan Likelihood Ratio Test.
Contoh penggunaan: Misalnya, teori menyatakan bahwa faktor x1 dan x2 mempengaruhi probabilitas seseorang membeli produk. Maka model logistik dibangun langsung dengan x1 dan x2, lalu diuji apakah kontribusi x2 benar-benar signifikan.
2. Exploratory (Pendekatan Eksploratori)
Pendekatan ini digunakan ketika peneliti belum memiliki teori yang pasti atau ingin mengeksplorasi hubungan potensial antar variabel.
Ciri-ciri:
• Model dibangun secara bertahap dengan tujuan menemukan kombinasi prediktor terbaik.
• Pemilihan variabel dilakukan berdasarkan kriteria statistik, seperti AIC, deviance, atau log-likelihood. Proses seleksi dilakukan melalui:
• Forward Selection: Mulai dari model kosong, satu per satu variabel dimasukkan jika signifkan.
• Backward Elimination: Mulai dari model penuh, variabel yang tidak signifkan dikeluarkan.
• Stepwise Selection: Gabungan dari keduanya, variabel dapat masuk dan keluar secara dinamis.
Tujuan:
Menemukan model yang parsimonious, yaitu cukup sederhana namun memiliki performa prediksi yang baik. Pemilihan antara pendekatan Confrmatory dan Exploratory bergantung pada tujuan penelitian. Jika ingin menguji hipotesis tertentu, gunakan pendekatan Confrmatory. Jika ingin menemukan model terbaik berdasarkan data, gunakan pendekatan Exploratory. Dalam praktiknya, kedua pendekatan ini sering digunakan secara komplementer: teori digunakan sebagai dasar, dan seleksi eksploratori dilakukan untuk menyempurnakan model.
Simulasi Data
Misalkan sebuah perusahaan e-commerce ingin memprediksi apakah pengunjung situs web akan membeli produk atau tidak berdasarkan karakteristik perilaku mereka saat mengakses situs. Adapun tujuan dari analisis ini adalah membangun model klasifikasi untuk memprediksi keputusan pembelian berdasarkan Waktu yang dihabiskan di halaman produk (dalam menit), Apakah pengguna berasal dari iklan media sosial (1 = ya, 0 = tidak), dan Skor ketertarikan produk yang dihitung dari klik & scroll (nilai normal).
Struktur variabelnya yaitu:
y: Faktor biner (1 = beli, 0 = tidak beli)
x1: Numerik (waktu dalam menit)
x2: Biner (asal dari iklan)
x3: Numerik (skor minat)
library(knitr)
library(dplyr)
library(ggplot2)
library(MASS)
## Warning: package 'MASS' was built under R version 4.4.3
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
library(caret)
## Warning: package 'caret' was built under R version 4.4.3
## Loading required package: lattice
library(pROC)
## Warning: package 'pROC' was built under R version 4.4.3
## Type 'citation("pROC")' for a citation.
##
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
##
## cov, smooth, var
library(DescTools)
## Warning: package 'DescTools' was built under R version 4.4.3
##
## Attaching package: 'DescTools'
## The following objects are masked from 'package:caret':
##
## MAE, RMSE
set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2) # waktu di halaman
x2 <- rbinom(n, 1, 0.4) # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1) # skor ketertarikan
# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3
# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))
# Simulasi variabel respon
y <- rbinom(n, 1, p)
# Buat data frame
df <- data.frame(
y = as.factor(y),
x1 = x1,
x2 = x2,
x3 = x3
)
head(df)
Pemilihan Model
Model Full
model_full <- glm(y ~ x1 + x2 + x3, data = df, family = binomial)
summary(model_full)
##
## Call:
## glm(formula = y ~ x1 + x2 + x3, family = binomial, data = df)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.2815 0.7045 -0.400 0.68947
## x1 0.5126 0.1597 3.210 0.00133 **
## x2 0.9536 0.5715 1.669 0.09518 .
## x3 0.7189 0.2530 2.842 0.00448 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 134.37 on 199 degrees of freedom
## Residual deviance: 112.44 on 196 degrees of freedom
## AIC: 120.44
##
## Number of Fisher Scoring iterations: 6
null_model <- glm(y ~ 1, data = df, family = binomial)
step_forward <- step(null_model, direction = "forward", scope = formula(model_full), trace = FALSE)
step_backward <- step(model_full, direction = "backward", trace = FALSE)
step_both <- step(null_model, direction = "both", scope = formula(model_full), trace = FALSE)
AIC(model_full, step_forward, step_backward, step_both)
pred_prob <- predict(step_both, type = "response")
roc_obj <- roc(df$y, pred_prob)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
plot(roc_obj, main = "Kurva ROC", col = "purple")
auc(roc_obj)
## Area under the curve: 0.7768
PseudoR2(step_both, which = c("CoxSnell", "Nagelkerke", "McFadden"))
## CoxSnell Nagelkerke McFadden
## 0.1038777 0.2123225 0.1632447
pred_class <- ifelse(pred_prob >= 0.5, 1, 0)
conf_matrix <- confusionMatrix(factor(pred_class), df$y, positive = "1")
conf_matrix
## Confusion Matrix and Statistics
##
## Reference
## Prediction 0 1
## 0 2 1
## 1 19 178
##
## Accuracy : 0.9
## 95% CI : (0.8498, 0.9378)
## No Information Rate : 0.895
## P-Value [Acc > NIR] : 0.4661224
##
## Kappa : 0.1442
##
## Mcnemar's Test P-Value : 0.0001439
##
## Sensitivity : 0.99441
## Specificity : 0.09524
## Pos Pred Value : 0.90355
## Neg Pred Value : 0.66667
## Prevalence : 0.89500
## Detection Rate : 0.89000
## Detection Prevalence : 0.98500
## Balanced Accuracy : 0.54483
##
## 'Positive' Class : 1
##
conf_matrix$byClass[c("Sensitivity", "Specificity")]
## Sensitivity Specificity
## 0.9944134 0.0952381
library(MASS)
library(broom)
library(DescTools)
Simulasi Data
set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2) # waktu di halaman
x2 <- rbinom(n, 1, 0.4) # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1) # skor ketertarikan
# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3
# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))
# Simulasi variabel respon
y <- rbinom(n, 1, p)
# Buat data frame
df <- data.frame(
y = as.factor(y),
x1 = x1,
x2 = x2,
x3 = x3
)
head(df)
Pembuatan Model
model1 <- glm(y ~ x1, data = df, family = binomial)
model2 <- glm(y ~ x1 + x2, data = df, family = binomial)
model3 <- glm(y ~ x1 + x2 + x3, data = df, family = binomial)
Perbandingan AIC dan Deviance
model_comp <- data.frame(
Model = c("Model 1", "Model 2", "Model 3"),
AIC = c(AIC(model1), AIC(model2), AIC(model3)),
Deviance = c(deviance(model1), deviance(model2), deviance(model3))
)
model_comp
anova(model1, model2, test = "LRT")
anova(model2, model3, test = "LRT")
Model yang kompleks sering memiliki AIC dan deviance yang lebih kecil, namun Model sederhana lebih mudah diinterpretasikan.
Jika penurunan AIC tidak signifikan, pilih model lebih sederhana.
Parsimony mencegah overfitting.
Rumus dan Penjelasan
Rumus AIC
\[ \text{AIC} = -2(\log L - k) = -2 \log L + 2k \]
Penjelasan:
AIC adalah ukuran untuk menilai model berdasarkan kombinasi antara
goodness-of-fit (melalui log-likelihood) dan kompleksitas
(melalui jumlah parameter \(k\)).
Semakin kecil AIC, semakin baik model tersebut secara keseluruhan karena
AIC menghukum model yang terlalu kompleks meskipun memiliki likelihood
tinggi.
Rumus Deviance
\[ D = -2 \left[\log L(\text{model}) - \log L(\text{model saturasi})\right] \]
Penjelasan: Deviance mengukur seberapa jauh model saat ini dibandingkan dengan model sempurna (model saturasi). Nilai deviance yang kecil menu
Rumus Likelihood-Ratio
\[ G^2 = -2 (\log L_0 - \log L_1) \]
Penjelasan:
Statistik Likelihood Ratio digunakan untuk menguji apakah penambahan
variabel dalam model secara signifikan meningkatkan kecocokan model.
Jika \(G^2\) besar dan p-value
kecil, maka model kompleks lebih baik dari model sederhana secara
statistik.
pred_prob <- predict(model3, type = "response")
pred_class <- factor(ifelse(pred_prob >= 0.5, 1, 0))
conf_matrix <- confusionMatrix(pred_class, df$y, positive = "1")
conf_matrix
## Confusion Matrix and Statistics
##
## Reference
## Prediction 0 1
## 0 2 1
## 1 19 178
##
## Accuracy : 0.9
## 95% CI : (0.8498, 0.9378)
## No Information Rate : 0.895
## P-Value [Acc > NIR] : 0.4661224
##
## Kappa : 0.1442
##
## Mcnemar's Test P-Value : 0.0001439
##
## Sensitivity : 0.99441
## Specificity : 0.09524
## Pos Pred Value : 0.90355
## Neg Pred Value : 0.66667
## Prevalence : 0.89500
## Detection Rate : 0.89000
## Detection Prevalence : 0.98500
## Balanced Accuracy : 0.54483
##
## 'Positive' Class : 1
##
Sensitivitas: Kemampuan model mendeteksi kelas positif secara benar (True Positive Rate)
\[ \text{Sensitivity} = \frac{\text{TP}}{\text{TP} + \text{FN}} \]
Spesifisitas: Kemampuan model mendeteksi kelas negatif secara benar (True Negative Rate)
\[ \text{Specificity} = \frac{\text{TN}}{\text{TN} + \text{FP}} \]
conf_matrix$byClass[c("Sensitivity", "Specificity")]
## Sensitivity Specificity
## 0.9944134 0.0952381
Kesimpulan
• Deviance yang kecil menunjukkan kecocokan model yang lebih baik.
• AIC yang rendah menunjukkan keseimbangan antara kecocokan dan kompleksitas.
• Likelihood Ratio Test mengevaluasi apakah model kompleks secara signifkan lebih baik.
• Tabel klasifkasi membantu menilai kinerja prediksi aktual vs prediksi model.
• Prinsip Parsimony mengutamakan model sederhana jika performanya mirip.
Kurva ROC adalah alat visual yang digunakan untuk mengevaluasi performa model klasifkasi biner. Kurva ini menunjukkan trade-of antara True Positive Rate (Sensitivity) dan False Positive Rate (1 - Specificity) pada berbagai threshold klasifkasi.
• Sumbu Y: Sensitivity = True Positive Rate = TP / (TP + FN)
• Sumbu X: 1 - Specifcity = False Positive Rate = FP / (FP + TN)
• Garis diagonal (dari kiri bawah ke kanan atas) menunjukkan performa acak (random guess).
• Kurva yang mendekati pojok kiri atas menunjukkan performa klasifkasi yang lebih baik.
Saat cut-off menurun, model mengklasifkasikan lebih banyak pengamatan sebagai positif:
• Sensitivitas naik
• Spesifisitas turun
Saat cut-off naik, model menjadi lebih konservatif:
• Sensitivitas turun
• Spesifisitas naik
Kurva ideal memiliki bentuk:
• Naik tajam secara vertikal hingga mencapai sensitivitas = 1
• Lalu bergerak secara horizontal menuju 1 - specifcity = 1
• Area under the curve (AUC) mendekati 1
• AUC = 0.5: model tidak lebih baik dari tebak acak
• AUC > 0.7: model cukup baik
• AUC > 0.9: model sangat baik
• AUC dikenal juga sebagai concordance index, yaitu probabilitas bahwa model memberikan nilai skor probabilitas yang lebih tinggi untuk kasus positif daripada kasus negatif.
• Untuk membandingkan performa beberapa model klasifkasi
• Untuk memilih threshold (cut-of) optimal berdasarkan kebutuhan aplikasi (misalnya: lebih penting menghindari false negative atau false positive)
library(pROC)
set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2) # waktu di halaman
x2 <- rbinom(n, 1, 0.4) # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1) # skor ketertarikan
# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3
# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))
# Simulasi variabel respon
y <- rbinom(n, 1, p)
data <- data.frame(y = as.factor(y), x1, x2, x3)
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
pred <- predict(model, type = "response")
roc_obj <- roc(data$y, pred)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
plot(roc_obj)
auc(roc_obj)
## Area under the curve: 0.7768
Untuk memilih threshold terbaik, kita bisa mengevaluasi sensitivitas dan spesifsitas pada berbagai cut-off.
thresholds <- seq(0.1, 0.9, by = 0.05)
results <- data.frame(Threshold = thresholds)
true_y <- factor(data$y, levels = c(0, 1))
results$Sensitivity <- sapply(thresholds, function(t) {
pred_class <- factor(ifelse(pred >= t, 1, 0), levels = c(0, 1))
cm <- table(Pred = pred_class, Obs = true_y)
TP <- ifelse("1" %in% rownames(cm) && "1" %in% colnames(cm), cm["1", "1"], 0)
FN <- ifelse("0" %in% rownames(cm) && "1" %in% colnames(cm), cm["0", "1"], 0)
if ((TP + FN) == 0) return(NA)
TP / (TP + FN)
})
results$Specificity <- sapply(thresholds, function(t) {
pred_class <- factor(ifelse(pred >= t, 1, 0), levels = c(0, 1))
cm <- table(Pred = pred_class, Obs = true_y)
TN <- ifelse("0" %in% rownames(cm) && "0" %in% colnames(cm), cm["0", "0"], 0)
FP <- ifelse("1" %in% rownames(cm) && "0" %in% colnames(cm), cm["1", "0"], 0)
if ((TN + FP) == 0) return(NA)
TN / (TN + FP)
})
print(results)
## Threshold Sensitivity Specificity
## 1 0.10 1.0000000 0.00000000
## 2 0.15 1.0000000 0.00000000
## 3 0.20 1.0000000 0.00000000
## 4 0.25 1.0000000 0.00000000
## 5 0.30 1.0000000 0.00000000
## 6 0.35 1.0000000 0.04761905
## 7 0.40 1.0000000 0.04761905
## 8 0.45 0.9944134 0.09523810
## 9 0.50 0.9944134 0.09523810
## 10 0.55 0.9888268 0.14285714
## 11 0.60 0.9832402 0.19047619
## 12 0.65 0.9832402 0.28571429
## 13 0.70 0.9664804 0.33333333
## 14 0.75 0.9553073 0.33333333
## 15 0.80 0.8938547 0.38095238
## 16 0.85 0.8212291 0.61904762
## 17 0.90 0.6480447 0.71428571
Cut-off optimal bisa dipilih berdasarkan:
Maksimum dari Sensitivity + Specificity
Atau mempertimbangkan trade-of sesuai tujuan aplikasi (misalnya: jika False Negative harus dihindari, maka prioritaskan sensitivitas tinggi)
ROC cocok saat proporsi kelas seimbang
Untuk data dengan kelas tidak seimbang, precision-recall curve bisa lebih informatif
Penjelasan Precision-Recall Curve
Kurva Precision-Recall (PR) adalah alat evaluasi performa model klasifikasi, khususnya sangat berguna saat bekerja dengan data yang tidak seimbang (class imbalance).
1. Definisi
Precision (Presisi): Proporsi prediksi positif yang benar-benar positif
\[ \text{Precision} = \frac{\text{TP}}{\text{TP} + \text{FP}} \]
Recall (Sensitivitas): Proporsi kasus positif yang berhasil diprediksi positif
\[ \text{Recall} = \frac{\text{TP}}{\text{TP} + \text{FN}} \]
2. Interpretasi
3. Area Under PR Curve
4. PR Curve vs ROC Curve
Aspek | ROC Curve | Precision-Recall Curve |
---|---|---|
Fokus | Semua kelas | Kelas positif saja |
Kuat di | Data seimbang | Data tidak seimbang |
Sumbu Y | Sensitivitas (Recall) | Precision |
Sumbu X | 1 - Spesifisitas | Recall |
5. Visualisasi PR Curve di R
library(PRROC)
## Warning: package 'PRROC' was built under R version 4.4.3
## Loading required package: rlang
set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2) # waktu di halaman
x2 <- rbinom(n, 1, 0.4) # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1) # skor ketertarikan
# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3
# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))
# Simulasi variabel respon
y <- rbinom(n, 1, p)
data <- data.frame(y = y, x1, x2, x3)
head(data)
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
prob <- predict(model, type = "response")
pr <- pr.curve(scores.class0 = prob[data$y == 1],
scores.class1 = prob[data$y == 0],
curve = TRUE)
plot(pr)
6. Catatan
• PR Curve sangat informatif untuk aplikasi seperti deteksi penipuan atau diagnosis penyakit langka.
• PR Curve digunakan saat:
– Kelas positif jauh lebih jarang daripada kelas negatif
– Tujuan aplikasi lebih mementingkan presisi terhadap kelas minoritas
Simulasi Data
set.seed(123)
n <- 200
x1 <- rnorm(n, mean = 5, sd = 2) # waktu di halaman
x2 <- rbinom(n, 1, 0.4) # asal dari iklan
x3 <- rnorm(n, mean = 0, sd = 1) # skor ketertarikan
# Membentuk prediktor linier
lin_pred <- -1 + 0.6 * x1 + 1.0 * x2 + 0.8 * x3
# Konversi ke probabilitas (fungsi logistik)
p <- 1 / (1 + exp(-lin_pred))
# Simulasi variabel respon
y <- rbinom(n, 1, p)
data <- data.frame(y = as.factor(y), x1, x2, x3)
head(data)
Model Logistik dan Null Model
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
model_null <- glm(y ~ 1, data = data, family = binomial)
Likelihood dan Rumus \[ R^2_{\text{Cox and Snell}} = 1 - \left( \frac{L_0}{L_M} \right)^{2/n} \]
\[ R^2_{\text{McFadden}} = 1 - \frac{\log L_M}{\log L_0} \]
Dengan:
Perhitungan Manual R-squared
logL0 <- logLik(model_null)
logLM <- logLik(model)
L0 <- exp(logL0)
LM <- exp(logLM)
n <- nobs(model)
cox_snell <- 1 - (L0 / LM)^(2 / n)
mcfadden <- 1 - (as.numeric(logLM) / as.numeric(logL0))
r2 <- data.frame(
R2_Cox_Snell = cox_snell,
R2_McFadden = mcfadden
)
r2
Perhitungan Otomatis dengan Package Tambahan
Menggunakan DescTools
if (!require(DescTools)) install.packages("DescTools"); library(DescTools)
PseudoR2(model, which = "all")
## McFadden McFaddenAdj CoxSnell Nagelkerke AldrichNelson
## 0.16324473 0.10370891 0.10387765 0.21232254 0.09883794
## VeallZimmermann Efron McKelveyZavoina Tjur AIC
## 0.24594790 0.14859709 0.31510399 0.14005131 120.43722088
## BIC logLik logLik0 G2
## 133.63049035 -56.21861044 -67.18644287 21.93566486
Interpretasi
Nilai McFadden R² sebesar 0.163 menunjukkan bahwa model memiliki kecocokan yang moderat.
Nilai Cox & Snell R² sebesar 0.1039 bersifat lebih konservatif dan memang tidak bisa mencapai 1 penuh karena keterbatasan metode ini. Nilai ini menunjukkan bahwa sekitar 10.4% variasi dalam data dijelaskan oleh model.
Nilai Nagelkerke R² sebesar 0.2123 merupakan versi koreksi dari Cox & Snell agar dapat mencapai maksimum 1. Dengan nilai > 0.2, ini memberikan indikasi bahwa model memiliki kecocokan yang layak
Distribusi multinomial adalah perluasan dari distribusi binomial untuk lebih dari dua kategori.
Jika \(X_1, X_2, \ldots, X_k\) menyatakan banyaknya kejadian dalam masing-masing dari \(k\) kategori, maka:
\[ P(X_1 = x_1, \ldots, X_k = x_k) = \frac{n!}{x_1! x_2! \ldots x_k!} p_1^{x_1} p_2^{x_2} \ldots p_k^{x_k} \]
dengan \(\sum_{i=1}^k x_i = n\) dan \(\sum_{i=1}^k p_i = 1\).
Sebuah survei dilakukan terhadap 10 orang yang diminta memilih satu dari tiga jenis buah favorit: Apel (A), Jeruk (J), dan Pisang (P).
Hasil survei:
Probabilitas teoretik preferensi:
Pertanyaannya: Berapa peluang bahwa dalam 10 orang akan ada 4 yang memilih Apel, 3 memilih Jeruk, dan 3 memilih Pisang?
Rumus Distribusi Multinomial
Distribusi peluang multinomial:
\[ P(X_1 = x_1, \ldots, X_k = x_k) = \frac{n!}{x_1! x_2! \ldots x_k!} p_1^{x_1} p_2^{x_2} \ldots p_k^{x_k} \]
Dengan:
Perhitungan Manual di R
n <- 10
x <- c(4, 3, 3)
p <- c(0.3, 0.4, 0.3)
# Hitung komponen-koefisien
faktorial_total <- factorial(n)
faktorial_x <- prod(factorial(x))
koefisien <- faktorial_total / faktorial_x
# Hitung peluang
peluang <- koefisien * prod(p^x)
peluang
## [1] 0.05878656
Probabilitas bahwa 4 orang memilih Apel, 3 Jeruk, dan 3 Pisang (dengan proporsi preferensi 0.3, 0.4, dan 0.3) adalah 0.058787. Distribusi multinomial digunakan untuk menghitung peluang dalam percobaan dengan beberapa kategori hasil. Rumus dasarnya merupakan generalisasi dari binomial untuk lebih dari dua kategori.
Model ini digunakan untuk memodelkan hubungan antara satu variabel respon kategorik (>2 kategori) dan satu atau lebih variabel prediktor.
Misalkan \(Y\) memiliki \(K\) kategori, dan kita pilih referensi (baseline) kategori \(K\), maka model logit untuk kategori \(j\) adalah:
\[ \log \left( \frac{P(Y = j)}{P(Y = K)} \right) = \beta_{0} + \beta_{1j}x_{1} + \cdots + \beta_{pj}x_{p} \]
untuk \(j = 1, 2, \ldots, K - 1\).
Baseline-category logit model adalah model regresi logistik untuk variabel respon kategorik dengan lebih dari dua kategori (nominal). Model ini menggunakan satu kategori sebagai acuan (baseline) dan membandingkan kategori lainnya terhadap baseline tersebut dalam bentuk logit:
\[ \log \left( \frac{\pi_j}{\pi_c} \right), \quad j = 1, \ldots, c - 1 \]
dengan:
Maka, terdapat sebanyak \((c - 1)\) fungsi logit.
Catatan: Kategori baseline bisa ditentukan secara eksplisit, tetapi default di R adalah kategori terakhir.
Model Regresi
Jika terdapat satu prediktor \(x\), maka bentuk umum model logit-nya adalah:
\[ \log \left( \frac{\pi_j}{\pi_c} \right) = \alpha_j + \beta_j x, \quad j = 1, \ldots, c - 1 \]
Contoh Kasus: 3 Kategori Respon
Misalkan respon \(Y\) memiliki tiga kategori: \(Y \in \{1, 2, 3\}\), dan kita gunakan kategori ke-3 sebagai baseline. Maka:
\[ \log \left( \frac{\pi_1}{\pi_3} \right) = \alpha_1 + \beta_1 x \]
\[ \log \left( \frac{\pi_2}{\pi_3} \right) = \alpha_2 + \beta_2 x \]
Terdapat dua model logit, satu untuk perbandingan kategori 1 dengan 3, dan satu lagi untuk kategori 2 dengan 3.
Relasi Antar Kategori
Jika ingin menghitung logit antara kategori 1 dan 2:
\[ \log \left( \frac{\pi_1}{\pi_2} \right) = \log \left( \frac{\pi_1 / \pi_3}{\pi_2 / \pi_3} \right) = \log \left( \frac{\pi_1}{\pi_3} \right) - \log \left( \frac{\pi_2}{\pi_3} \right) \] \[ = (\alpha_1 + \beta_1 x) - (\alpha_2 + \beta_2 x) = (\alpha_1 - \alpha_2) + (\beta_1 - \beta_2)x \] Model Baseline-category Logit:
multinom()
dari
package nnet, dan kategori baseline bisa ditentukan
dengan relevel()
.Estimasi dilakukan menggunakan metode maximum likelihood dengan algoritma iteratif seperti Newton-Raphson.
Log-likelihood:
\[ l(\beta) = \sum_{i=1}^n \sum_{j=1}^K y_{ij} \log (\pi_{ij}) \]
dengan \(\pi_{ij} = P(Y = j | x_i)\) dan \(y_{ij} = 1\) jika \(Y_i = j\), 0 jika tidak.
Sebuah institusi pendidikan ingin memahami faktor-faktor yang memengaruhi preferensi gaya belajar siswa, yaitu Visual, Auditory, atau Kinesthetic. Sebanyak 150 siswa disurvei, dan dikumpulkan data sebagai berikut:
LearningStyle: Gaya belajar yang dipilih (Visual, Auditory, Kinesthetic)
GradeLevel: Tingkatan kelas siswa (Junior, Senior, Freshman)
StudyHours: Rata-rata jam belajar mandiri per minggu
Extracurriculars: Jumlah kegiatan ekstrakurikuler yang diikuti
Tujuan analisis ini adalah untuk mengetahui bagaimana GradeLevel, StudyHours, dan Extracurriculars memengaruhi preferensi gaya belajar siswa.
set.seed(123)
n <- 150
GradeLevel <- sample(c("Freshman", "Junior", "Senior"), n, replace = TRUE)
StudyHours <- round(rnorm(n, mean = 10, sd = 3))
Extracurriculars <- sample(0:5, n, replace = TRUE)
# Simulasikan gaya belajar berdasarkan grade level
LearningStyle <- sapply(GradeLevel, function(level) {
if (level == "Freshman") {
sample(c("Visual", "Auditory", "Kinesthetic"), size = 1, prob = c(0.5, 0.3, 0.2))
} else if (level == "Junior") {
sample(c("Visual", "Auditory", "Kinesthetic"), size = 1, prob = c(0.3, 0.4, 0.3))
} else {
sample(c("Visual", "Auditory", "Kinesthetic"), size = 1, prob = c(0.2, 0.3, 0.5))
}
})
# Buat data frame
df <- data.frame(
LearningStyle = factor(LearningStyle),
GradeLevel = factor(GradeLevel),
StudyHours,
Extracurriculars
)
# Tetapkan baseline untuk analisis (misal: Visual sebagai baseline)
df$LearningStyle <- relevel(df$LearningStyle, ref = "Visual")
head(df)
library(nnet)
model_mnlogit <- multinom(LearningStyle ~ GradeLevel + StudyHours + Extracurriculars, data = df)
## # weights: 18 (10 variable)
## initial value 164.791843
## iter 10 value 150.818303
## final value 150.779592
## converged
summary(model_mnlogit)
## Call:
## multinom(formula = LearningStyle ~ GradeLevel + StudyHours +
## Extracurriculars, data = df)
##
## Coefficients:
## (Intercept) GradeLevelJunior GradeLevelSenior StudyHours
## Auditory -1.313043 1.176353 0.9339879 0.09852580
## Kinesthetic -1.977086 1.818910 2.1623436 0.03622559
## Extracurriculars
## Auditory 0.0101990
## Kinesthetic 0.2123591
##
## Std. Errors:
## (Intercept) GradeLevelJunior GradeLevelSenior StudyHours
## Auditory 0.8645914 0.5164149 0.5376260 0.07033911
## Kinesthetic 0.9351475 0.5908354 0.5904723 0.07194458
## Extracurriculars
## Auditory 0.1330148
## Kinesthetic 0.1377687
##
## Residual Deviance: 301.5592
## AIC: 321.5592
z <- summary(model_mnlogit)$coefficients / summary(model_mnlogit)$standard.errors
pval <- 2 * (1 - pnorm(abs(z)))
round(pval, 4)
## (Intercept) GradeLevelJunior GradeLevelSenior StudyHours
## Auditory 0.1288 0.0227 0.0823 0.1613
## Kinesthetic 0.0345 0.0021 0.0003 0.6146
## Extracurriculars
## Auditory 0.9389
## Kinesthetic 0.1232
Interpretasi:
• Koefisien untuk kategori “Auditory” dan “Kinesthetic” dibandingkan dengan baseline “Visual”
• Nilai p-value kecil (<0.05) menunjukkan variabel tersebut signifkan memengaruhi preferensi gaya belajar.
df$Predicted <- predict(model_mnlogit)
table(Predicted = df$Predicted, Actual = df$LearningStyle)
## Actual
## Predicted Visual Auditory Kinesthetic
## Visual 18 12 4
## Auditory 10 26 17
## Kinesthetic 11 18 34
Model regresi logistik multinomial berhasil digunakan untuk:
• Menganalisis hubungan antara atribut sisawa dan preferensi gaya belajar
• Mengetahui faktor signifkan yang memengaruhi pilihan
• Memungkinkan prediksi gaya belajar yang dipilih oleh siswa berdasarkan karakteristiknya
Kita gunakan dataset iris untuk memodelkan jenis spesies bunga berdasarkan lebar dan panjang kelopak.
data(iris)
iris <- iris %>% mutate(Species = relevel(Species, ref = "setosa"))
model <- multinom(Species ~ Petal.Length + Petal.Width, data = iris)
## # weights: 12 (6 variable)
## initial value 164.791843
## iter 10 value 12.657828
## iter 20 value 10.374056
## iter 30 value 10.330881
## iter 40 value 10.306926
## iter 50 value 10.300057
## iter 60 value 10.296452
## iter 70 value 10.294046
## iter 80 value 10.292029
## iter 90 value 10.291154
## iter 100 value 10.289505
## final value 10.289505
## stopped after 100 iterations
summary(model)
## Call:
## multinom(formula = Species ~ Petal.Length + Petal.Width, data = iris)
##
## Coefficients:
## (Intercept) Petal.Length Petal.Width
## versicolor -22.79944 6.92122 7.878496
## virginica -67.82521 12.64721 18.261016
##
## Std. Errors:
## (Intercept) Petal.Length Petal.Width
## versicolor 44.3859 37.58715 81.00888
## virginica 46.3939 37.65702 81.09482
##
## Residual Deviance: 20.57901
## AIC: 32.57901
Interpretasi Koefisien
z <- summary(model)$coefficients / summary(model)$standard.errors
p_values <- 2 * (1 - pnorm(abs(z)))
round(p_values, 4)
## (Intercept) Petal.Length Petal.Width
## versicolor 0.6075 0.8539 0.9225
## virginica 0.1438 0.7370 0.8218
Nilai p-value menunjukkan apakah variabel prediktor berpengaruh signifkan terhadap log odds dibanding baseline category.
Prediksi dan Visualisasi
iris$predicted <- predict(model, newdata = iris)
table(Predicted = iris$predicted, Actual = iris$Species)
## Actual
## Predicted setosa versicolor virginica
## setosa 50 0 0
## versicolor 0 47 3
## virginica 0 3 47
ggplot(iris, aes(x = Petal.Length, y = Petal.Width, color = predicted)) +
geom_point(size = 2) +
labs(title = "Multinomial Logistic Regression Predictions",
x = "Petal Length", y = "Petal Width") +
theme_minimal()
Kesimpulan Multinomial logistic regression efektif digunakan untuk klasifikasi kategori nominal. Model ini memberikan estimasi log-odds terhadap baseline dan prediksi kategori baru berdasarkan kombinasi prediktor.
Regresi logistik ordinal digunakan ketika variabel respon \(Y\) bersifat ordinal (memiliki urutan), misalnya tingkat kepuasan: Rendah, Sedang, Tinggi.
Model ini berbeda dengan:
• Regresi logistik biner: hanya 2 kategori
• Regresi logistik multinomial: kategori > 2 tetapi tidak berurutan
Model yang digunakan adalah Cumulative Logit Model dengan asumsi proportional odds:
\[ \log \left( \frac{P(Y \leq j)}{P(Y > j)} \right) = \alpha_j + \beta x \]
Untuk \(c\) kategori, terdapat \((c - 1)\) model logit kumulatif.
Koefisien \(\beta\) menjelaskan efek \(x\) terhadap kemungkinan berada pada kategori yang lebih rendah atau sama.
Odds ratio ditulis:
\[ \text{OR} = e^\beta \]
Misalkan seorang dosen ingin mengetahui apakah semakin tinggi intensitas belajar (jam belajar per minggu) mahasiswa berhubungan dengan tingkat pemahaman mereka terhadap materi.
Variabel yang digunakan:
1: Rendah
2: Sedang
3: Tinggi
set.seed(123)
n <- 200
study_hours <- round(runif(n, 1, 20)) # Jam belajar antara 1 hingga 20 jam
understanding <- cut(3 + 0.3*study_hours + rnorm(n),
breaks = c(-Inf, 6.5, 9.5, Inf),
labels = c("Rendah", "Sedang", "Tinggi"),
ordered_result = TRUE)
df <- data.frame(understanding, study_hours)
head(df)
model_ord <- polr(understanding ~ study_hours, data = df, Hess = TRUE)
summary(model_ord)
## Call:
## polr(formula = understanding ~ study_hours, data = df, Hess = TRUE)
##
## Coefficients:
## Value Std. Error t value
## study_hours 0.5241 0.06772 7.74
##
## Intercepts:
## Value Std. Error t value
## Rendah|Sedang 6.1901 0.8153 7.5923
## Sedang|Tinggi 11.9717 1.3113 9.1299
##
## Residual Deviance: 162.6821
## AIC: 168.6821
(ctable <- coef(summary(model_ord)))
## Value Std. Error t value
## study_hours 0.5240819 0.06771506 7.739519
## Rendah|Sedang 6.1900621 0.81531258 7.592256
## Sedang|Tinggi 11.9717293 1.31126286 9.129923
p <- pnorm(abs(ctable[, "t value"]), lower.tail = FALSE) * 2
(ctable <- cbind(ctable, "p value" = round(p, 4)))
## Value Std. Error t value p value
## study_hours 0.5240819 0.06771506 7.739519 0
## Rendah|Sedang 6.1900621 0.81531258 7.592256 0
## Sedang|Tinggi 11.9717293 1.31126286 9.129923 0
newdata <- data.frame(study_hours = 5:9)
predict(model_ord, newdata = newdata, type = "probs")
## Rendah Sedang Tinggi
## 1 0.9726059 0.02730723 0.0000868432
## 2 0.9545909 0.04526244 0.0001466614
## 3 0.9256346 0.07411775 0.0002476725
## 4 0.8805245 0.11905726 0.0004182246
## 5 0.8135622 0.18573168 0.0007061392
Model cumulative logit mengasumsikan efek prediktor sama untuk setiap cutoff. Jika tidak, pertimbangkan model non-proportional odds seperti generalized ordinal model.
Selain cumulative logit, model ordinal lainnya:
• Adjacent-category logit
• Continuation-ratio (sequential) logit
Model tersebut dapat digunakan saat asumsi proportional odds tidak terpenuhi
polr()
dari package
MASS.Jika diperlukan validasi lanjut, bisa digunakan uji devian atau likelihood ratio test.
Model regresi logistik ordinal yang paling umum digunakan adalah Cumulative Logit Model dengan asumsi Proportional Odds.
Asumsi ini dikenal juga sebagai asumsi paralelisme (parallel lines assumption).
Asumsi paralelisme menyatakan bahwa koefisien regresi (\(\beta\)) sama untuk setiap kategori kumulatif dari variabel respon.
Bentuk umum model:
\[ \log\left(\frac{P(Y \leq j)}{P(Y > j)}\right) = \alpha_j + \beta x \]
Untuk \(j = 1, 2, \dots, c - 1\):
Visualisasi: Dalam asumsi paralelisme, kurva logit kumulatif dari tiap kategori terhadap prediktor akan memiliki kemiringan yang sama (paralel), hanya berbeda posisi (intercept).
Konsekuensi Pelanggaran Asumsi
Jika asumsi ini tidak terpenuhi:
Pengujian Asumsi Paralelisme
Untuk memeriksa validitas asumsi, dapat digunakan:
Kesimpulan
• Asumsi paralelisme penting untuk validitas model cumulative logit.
• Menyederhanakan interpretasi karena efek prediktor konstan.
• Jika tidak terpenuhi, gunakan model ordinal alternatif.
Analisis data kategorikal memegang peranan penting dalam statistika terapan karena banyak kejadian di dunia nyata menghasilkan data dalam kategori, seperti jenis kelamin, status pekerjaan, tingkat pendidikan, preferensi konsumen, atau diagnosis medis. Data yang berbentuk kategori ini biasanya dianalisis dengan menggunakan tabel kontingensi, model log-linier, dan model regresi logistik. Setiap metode memiliki keuntungan dan kerugian yang bervariasi tergantung pada tujuan analisis serta struktur data.
Tabel kontingensi berfungsi sebagai langkah awal untuk mengeksplorasi hubungan antara dua atau lebih variabel kategorikal. Contohnya, dalam penelitian mengenai dampak obat terhadap serangan jantung, tabel kontingensi dapat menunjukkan jumlah pasien yang mengalami atau tidak mengalami serangan jantung berdasarkan obat yang digunakan. Tabel ini berperan dalam mengidentifikasi pola-pola awal serta menghitung ukuran asosiasi seperti odds ratio, risk ratio, dan statistik chi-square untuk menguji independensi antara variabel.
Namun, untuk menciptakan model statistik yang dapat mengontrol efek dari berbagai variabel dan interaksinya secara bersamaan, model log-linier menjadi sangat bermanfaat. Model log-linier merupakan tipe khusus dari Generalized Linear Model (GLM) yang diterapkan pada frekuensi sel dalam tabel kontingensi dan mengasumsikan distribusi Poisson. Berbeda dengan regresi logistik, model log-linier tidak menentukan variabel mana yang dependen dan mana yang independen, melainkan memperlakukannya secara simetris. Model ini lebih sesuai bila tujuan analisis adalah untuk memahami struktur asosiasi atau independensi antar variabel, alih-alih untuk tujuan prediksi.
Struktur model log-linier ditentukan oleh efek utama dari masing-masing variabel serta interaksi di antara mereka. Sebagai contoh, dalam tabel kontingensi tiga arah (seperti: jenis kelamin, status merokok, dan penyakit paru-paru), model ini mampu mengidentifikasi apakah interaksi antara dua variabel cukup untuk menjelaskan data, atau apakah interaksi tiga arah diperlukan untuk memahami struktur asosiasi yang ada. Penyesuaian model bisa dilakukan melalui metode uji rasio kemungkinan untuk membandingkan model yang lebih sederhana dengan yang lebih rumit.
Sebaliknya, regresi logistik adalah metode yang paling umum dipakai ketika satu variabel kategorikal secara spesifik dianggap sebagai variabel dependen (misalnya, kejadian penyakit: ya/tidak) dan satu atau lebih variabel kategorik maupun numerik berfungsi sebagai prediktor. Model ini mengkaji logit dari probabilitas kejadian (atau log odds), dan sangat berguna dalam penelitian observasional serta eksperimental untuk menjelaskan atau meramalkan kemungkinan suatu hasil. Regresi logistik juga memiliki variasi untuk outcome kategorik yang lebih dari dua kelas, termasuk regresi logistik multinomial dan regresi logistik ordinal.
Oleh karena itu, meskipun ketiga metode ini beroperasi pada data kategorikal, tabel kontingensi bersifat deskriptif, model log-linier bersifat eksploratif dalam hubungan simetris, sedangkan regresi logistik adalah alat prediktif untuk outcome kategorikal. Pemilihan metode yang tepat sangat bergantung pada apakah analisis berfokus pada deskripsi, eksplorasi struktur, atau prediksi hasil berdasarkan variabel penjelas. Paduan ketiga metode ini sering dijumpai dalam praktik untuk memperoleh pemahaman yang lebih mendalam tentang data kategorikal yang dianalisis.
Ringkasan Dalam analisis data kategorikal, terdapat beberapa metode statistik yang umum dipakai, di antaranya: 1. Tabel Kontingensi: menampilkan frekuensi gabungan dari dua atau lebih variabel kategorikal. 2. Model Log-linier: digunakan untuk memodelkan struktur asosiasi di dalam tabel kontingensi tanpa menganggap adanya variabel dependen. 3. Model Regresi Logistik: digunakan untuk memodelkan probabilitas dari kategori variabel dependen berdasarkan variabel independen.
Meskipun ketiga metode ini dapat diaplikasikan pada data kategorikal, pendekatan dan cara interpretasinya sangatlah berbeda.
Tabel Kontingensi
Tabel kontingensi menyajikan jumlah frekuensi dari kombinasi kategori antar variabel.
Contoh tabel 2x2:
# Membuat matriks dari data
table_data <- matrix(c(28, 656, 18, 658),
nrow = 2,
byrow = TRUE,
dimnames = list(Obat = c("Placebo", "Aspirin"),
Infeksi = c("Ya", "Tidak")))
table_data
## Infeksi
## Obat Ya Tidak
## Placebo 28 656
## Aspirin 18 658
Tabel kontingensi bersifat deskriptif dan tidak melibatkan pemodelan probabilitas
Model Loglinear
Model loglinear memodelkan logaritma dari ekspektasi frekuensi sel dalam tabel kontingensi.
\[ \log(\mu_{ij}) = \mu + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]
Cocok untuk menyelidiki asosiasi dan independensi antar variabel.
Tidak membedakan antara variabel respon dan penjelas.
Umumnya digunakan pada tabel multi-dimensi (2 arah, 3 arah, dst).
library(MASS)
loglm(~ Obat * Infeksi, data = table_data)
## Call:
## loglm(formula = ~Obat * Infeksi, data = table_data)
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0 0 1
## Pearson 0 0 1
Model Regresi Logistik
Model regresi logistik biner:
\[ \log\left( \frac{p}{1 - p} \right) = \beta_0 + \beta_1 x \]
Contoh:
data_glm <- data.frame(
Infeksi = c(1, 0, 1, 0),
Obat = factor(c("Placebo", "Placebo", "Aspirin", "Aspirin")),
Frek = c(28, 656, 18, 658)
)
head(data_glm)
model_logit <- glm(Infeksi ~ Obat, weights = Frek, family = binomial, data = data_glm)
summary(model_logit)
##
## Call:
## glm(formula = Infeksi ~ Obat, family = binomial, data = data_glm,
## weights = Frek)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.5988 0.2389 -15.065 <2e-16 ***
## ObatPlacebo 0.4449 0.3071 1.449 0.147
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 401.99 on 3 degrees of freedom
## Residual deviance: 399.85 on 2 degrees of freedom
## AIC: 403.85
##
## Number of Fisher Scoring iterations: 6
Perbandingan Ketiga Pendekatan
Aspek | Tabel Kontingensi | Model Loglinear | Regresi Logistik |
---|---|---|---|
Tujuan | Deskripsi frekuensi | Deteksi asosiasi | Prediksi probabilitas |
Variabel dependen | Tidak ada | Tidak ada (intersit) | Ada (eksplisit) |
Distribusi | Tidak diasumsikan | Poisson (frekuensi sel) | Binomial (probabilitas) |
Bentuk Model | Tidak ada | GLM: log(μ) ~ efek | GLM: logit(p) ~ prediktor |
Cocok untuk | Eksplorasi awal | Tabel > 2 variabel | Studi prediktif |
Tabel kontingensi menyajikan frekuensi dari kombinasi kategori antar dua atau lebih variabel. Misal:
# Contoh tabel 2x2
table_data <- matrix(c(28, 656, 18, 658),
nrow = 2,
byrow = TRUE,
dimnames = list(Obat = c("Placebo", "Aspirin"),
Infeksi = c("Ya", "Tidak")))
table_data
## Infeksi
## Obat Ya Tidak
## Placebo 28 656
## Aspirin 18 658
Model log-linier untuk tabel I x J dapat dituliskan: \[ \log(\mu_{ij}) = \mu + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]
Model saturated atau model penuh menyertakan seluruh efek utama dan interaksi:
• Cocok sempurna terhadap data
• Tidak mengasumsikan independensi antar variabel
Contoh formulasi untuk tabel 2x2:
# Data
library(MASS)
data <- matrix(c(28, 656, 18, 658),
nrow = 2,
byrow = TRUE,
dimnames = list(Obat = c("Placebo", "Aspirin"),
Infeksi = c("Ya", "Tidak")))
ftable(data)
## Infeksi Ya Tidak
## Obat
## Placebo 28 656
## Aspirin 18 658
Model saturated dapat dipasang dengan loglm
dari package
{MASS}
:
model_saturated <- loglm(~ Obat * Infeksi, data = data)
summary(model_saturated)
## Formula:
## ~Obat * Infeksi
## attr(,"variables")
## list(Obat, Infeksi)
## attr(,"factors")
## Obat Infeksi Obat:Infeksi
## Obat 1 0 1
## Infeksi 0 1 1
## attr(,"term.labels")
## [1] "Obat" "Infeksi" "Obat:Infeksi"
## attr(,"order")
## [1] 1 1 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0 0 1
## Pearson 0 0 1
Model independen mengasumsikan bahwa tidak ada interaksi antara variabel:
\[ \log(\mu_{ij}) = \mu + \lambda^T_i + \lambda^B_j \] Model ini menguji hipotesis bahwa variabel X dan Y saling independen.
model_indep <- loglm(~ Obat + Infeksi, data = data)
summary(model_indep)
## Formula:
## ~Obat + Infeksi
## attr(,"variables")
## list(Obat, Infeksi)
## attr(,"factors")
## Obat Infeksi
## Obat 1 0
## Infeksi 0 1
## attr(,"term.labels")
## [1] "Obat" "Infeksi"
## attr(,"order")
## [1] 1 1
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 2.147353 1 0.1428159
## Pearson 2.129972 1 0.1444434
Odds ratio untuk tabel 2x2:
\[ OR = \frac{n_{11}n_{22}}{n_{12}n_{21}} \]
Interpretasi nilai OR:
Dalam model saturated:
# Estimasi odds ratio dan log-odds
logOR <- log((data[1,1] * data[2,2]) / (data[1,2] * data[2,1]))
logOR
## [1] 0.4448769
Perbandingan antar model dilakukan dengan menggunakan statistik deviance (G²) atau likelihood ratio test.
anova(model_indep, model_saturated)
## LR tests for hierarchical log-linear models
##
## Model 1:
## ~Obat + Infeksi
## Model 2:
## ~Obat * Infeksi
##
## Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1 2.147353 1
## Model 2 0.000000 0 2.147353 1 0.14282
## Saturated 0.000000 0 0.000000 0 1.00000
data_survey <- matrix(
c(25, 25, 12, 0, 1, 3),
nrow = 2,
byrow = TRUE,
dimnames = list(
Perilaku = c("Dikontrol", "Infeksi"),
Rokok = c("0", "1-24", ">25")
)
)
ftable(data_survey)
## Rokok 0 1-24 >25
## Perilaku
## Dikontrol 25 25 12
## Infeksi 0 1 3
loglm(~ Perilaku + Rokok, data = data_survey)
## Call:
## loglm(formula = ~Perilaku + Rokok, data = data_survey)
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 6.690106 2 0.03525835
## Pearson 6.956203 2 0.03086595
Model log-linear adalah model yang digunakan untuk menganalisis hubungan antara dua atau lebih variabel kategorik yang disajikan dalam tabel kontingensi. Model ini mengasumsikan bahwa logaritma dari nilai ekspektasi frekuensi sel (\(\mu_{ij}\)) dapat dinyatakan sebagai penjumlahan efek variabel dan (bila perlu) interaksinya. Untuk tabel 2x2
\[
\log(\mu_{ij}) = \lambda + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij}
\]
Sistem Persamaan Model Log-Linear
\[ \begin{align*} \log(\mu_{11}) &= \lambda + \lambda^A_1 + \lambda^B_1 + \lambda^{AB}_{11} \\ \log(\mu_{12}) &= \lambda + \lambda^A_1 + \lambda^B_2 + \lambda^{AB}_{12} \\ \log(\mu_{21}) &= \lambda + \lambda^A_2 + \lambda^B_1 + \lambda^{AB}_{21} \\ \log(\mu_{22}) &= \lambda + \lambda^A_2 + \lambda^B_2 + \lambda^{AB}_{22} \end{align*} \]
Constraint Sum-to-Zero
\[ \begin{aligned} &\lambda^A_1 + \lambda^A_2 = 0 \\ &\lambda^B_1 + \lambda^B_2 = 0 \\ &\lambda^{AB}_{11} + \lambda^{AB}_{12} + \lambda^{AB}_{21} + \lambda^{AB}_{22} = 0 \end{aligned} \]
Rumus Estimasi Parameter dengan Sum-to-Zero Constraint
\[ \lambda^A_1 = \frac{1}{2} \left[ (\log \mu_{11} + \log \mu_{12}) - (\log \mu_{21} + \log \mu_{22}) \right] \]
\[ \lambda^B_1 = \frac{1}{2} \left[ (\log \mu_{11} + \log \mu_{21}) - (\log \mu_{12} + \log \mu_{22}) \right] \]
\[ \lambda^{AB}_{12} = \frac{1}{4} \left[ \log \mu_{12} - \log \mu_{11} - \log \mu_{22} + \log \mu_{21} \right] \]
Diberikan data:
Infeksi | Tidak Infeksi | |
---|---|---|
Placebo | 28 | 656 |
Aspirin | 18 | 658 |
Model log-linear pada tabel 2x2:
\[ \log(\mu_{ij}) = \lambda + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]
dengan constraint sum-to-zero:
\[ \sum_i \lambda^A_i = 0,\quad \sum_j \lambda^B_j = 0,\quad \sum_{i,j} \lambda^{AB}_{ij} = 0 \]
Misalkan:
Observasi:
\[ \begin{aligned} \log(\mu_{11}) &= \lambda + \lambda^A_1 + \lambda^B_1 + \lambda^{AB}_{11} \\ \log(\mu_{12}) &= \lambda + \lambda^A_1 + \lambda^B_2 + \lambda^{AB}_{12} \\ \log(\mu_{21}) &= \lambda + \lambda^A_2 + \lambda^B_1 + \lambda^{AB}_{21} \\ \log(\mu_{22}) &= \lambda + \lambda^A_2 + \lambda^B_2 + \lambda^{AB}_{22} \end{aligned} \]
Dengan constraint sum-to-zero:
\[ \begin{aligned} \lambda^A_1 + \lambda^A_2 &= 0, \\ \lambda^B_1 + \lambda^B_2 &= 0, \\ \lambda^{AB}_{11} + \lambda^{AB}_{12} + \lambda^{AB}_{21} + \lambda^{AB}_{22} &= 0 \end{aligned} \]
Langkah-langkah:
1. Hitung rata-rata log frekuensi sel:
\[ \lambda = \frac{1}{4} \sum_{i=1}^{2} \sum_{j=1}^{2} \log(n_{ij}) \]
\[ = \frac{1}{4} (\log(28) + \log(656) + \log(18) + \log(658)) \]
\[ = \frac{1}{4} (3.3322 + 6.4850 + 2.8904 + 6.4887) \]
\[ = 4.7991 \]
2. Efek utama A (Jenis Obat):
\[ \lambda^A_1 = \frac{1}{2} [(\log(28) + \log(656)) - (\log(18) + \log(658))] \]
\[ = \frac{1}{2} [(3.3322 + 6.4850) - (2.8904 + 6.4887)] \]
\[ = \frac{1}{2} (0.4381) = 0.2190 \]
\[ \lambda^A_2 = -0.2190 \]
3. Efek utama B (Kategori Infeksi):
\[ \lambda^B_1 = \frac{1}{2} [(\log(28) + \log(18)) - (\log(656) + \log(658))] \]
\[ = \frac{1}{2} [(3.3322 + 2.8904) - (6.4850 + 6.4887)] \]
\[ = \frac{1}{2} (6.2226 - 12.9737) = \frac{1}{2} (-6.7511) = -3.3756 \]
\[ \lambda^B_2 = +3.3756 \]
4. Efek interaksi:
\[ \lambda^{AB}_{11} = \frac{1}{4} [\log(28) - \log(656) - \log(18) + \log(658)] \]
\[ = \frac{1}{4} [3.3322 - 6.4850 - 2.8904 + 6.4887] \]
\[ = \frac{1}{4} (0.4455) = 0.1114 \]
\[ \lambda^{AB}_{12} = -0.1114 \quad \lambda^{AB}_{21} = -0.1114 \quad \lambda^{AB}_{22} = +0.1114 \]
Ringkasan parameter:
\[ \lambda = 4.7991 \]
\[ \lambda^A_1 = 0.2190, \quad \lambda^A_2 = -0.2190 \]
\[ \lambda^B_1 = -3.3756, \quad \lambda^B_2 = 3.3756 \]
\[
\lambda^{AB}_{11} = 0.1114, \quad
\lambda^{AB}_{12} = -0.1114, \quad
\lambda^{AB}_{21} = -0.1114, \quad
\lambda^{AB}_{22} = 0.1114
\]
## 15.6 Hitung Odds Ratio dan Interval Kepercayaan
\[ \text{OR} = \frac{n_{11}n_{22}}{n_{12}n_{21}} = \frac{28 \times 658}{656 \times 18} = \frac{18424}{11808} = 1.5609 \]
Log odds ratio:
\[ \log(\text{OR}) = \log(1.5609) = 0.4455 \]
Standard error (SE):
\[ SE = \sqrt{\frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}}} \]
\[ SE = \sqrt{ \frac{1}{28} + \frac{1}{656} + \frac{1}{18} + \frac{1}{658} } = \sqrt{0.0357 + 0.0015 + 0.0556 + 0.0015} = \sqrt{0.0943} \] \[ = 0.3071 \]
95% Confidence Interval for log(OR):
\[ \log(OR) \pm 1.96 \times SE = 0.4455 \pm 1.96 \times 0.3071 \]
\[ = (0.4455 - 0.6019,\; 0.4455 + 0.6019) = (-0.1564,\; 1.0474) \]
Back-transform to get CI for OR:
\[ \text{Lower} = \exp(-0.1564) = 0.8555 \]
\[ \text{Upper} = \exp(1.0474) = 2.8501 \]
Jadi, OR = 1.56 (95% CI: 0.86 – 2.85)
# Data 2x2
tabel <- matrix(c(28, 656, 18, 658), nrow = 2, byrow = TRUE)
colnames(tabel) <- c("Infeksi", "Tidak")
rownames(tabel) <- c("Placebo", "Aspirin")
tabel
## Infeksi Tidak
## Placebo 28 656
## Aspirin 18 658
# Konversi jadi data.frame untuk model GLM
data <- as.data.frame(as.table(tabel))
colnames(data) <- c("Obat", "Infeksi", "Freq")
data
# Fit model tanpa interaksi
fit_no_inter <- glm(Freq ~ Obat + Infeksi, family = poisson(), data = data)
summary(fit_no_inter)
##
## Call:
## glm(formula = Freq ~ Obat + Infeksi, family = poisson(), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.14136 0.14989 20.958 <2e-16 ***
## ObatAspirin -0.01176 0.05423 -0.217 0.828
## InfeksiTidak 3.35219 0.15000 22.348 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 1485.5614 on 3 degrees of freedom
## Residual deviance: 2.1474 on 1 degrees of freedom
## AIC: 34.713
##
## Number of Fisher Scoring iterations: 4
# Fit model dengan interaksi
fit_inter <- glm(Freq ~ Obat * Infeksi, family = poisson(), data = data)
summary(fit_inter)
##
## Call:
## glm(formula = Freq ~ Obat * Infeksi, family = poisson(), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.3322 0.1890 17.632 <2e-16 ***
## ObatAspirin -0.4418 0.3021 -1.462 0.144
## InfeksiTidak 3.1540 0.1930 16.344 <2e-16 ***
## ObatAspirin:InfeksiTidak 0.4449 0.3071 1.449 0.147
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 1.4856e+03 on 3 degrees of freedom
## Residual deviance: 1.3323e-15 on 0 degrees of freedom
## AIC: 34.565
##
## Number of Fisher Scoring iterations: 3
Parameter utama (intercept) menunjukkan rata-rata log frekuensi sel referensi (Placebo + Infeksi).
Efek “Obat” dan “Infeksi” menunjukkan perbedaan log-frekuensi antara kategori. Efek infeksi signifikan; efek obat tidak.
Efek interaksi tidak signifikan, namun positif, menunjukkan bahwa efek Aspirin terhadap infeksi mungkin berbeda tergantung status infeksi, meski bukti statistiknya belum kuat
Sebuah studi meneliti hubungan antara kejadian Myocardial Infarction (MI) dan tingkat merokok (Smoking Level):
0 batang | 1–24 batang | >25 batang | |
---|---|---|---|
kontrol | 25 | 25 | 12 |
Infeksi Myocardial | 0 | 1 | 3 |
Bentuk umum model log-linear untuk tabel 2x3 (dengan sum-to-zero constraint):
\[ \log(\mu_{ij}) = \lambda + \lambda^A_i + \lambda^B_j + \lambda^{AB}_{ij} \]
dengan:
Constraint: \[ \sum_i \lambda^A_i = 0,\quad \sum_j \lambda^B_j = 0,\quad \sum_i \lambda^{AB}_{ij} = 0,\quad \sum_j \lambda^{AB}_{ij} = 0 \]
Secara eksplisit:
\[ \log(\mu_{ij}) = \lambda + \begin{cases} \lambda^A_1 & \text{(Kontrol)} \\ \lambda^A_2 & \text{(Infeksi Myocardial)} \end{cases} + \begin{cases} \lambda^B_1 & \text{(Merokok = 0)} \\ \lambda^B_2 & \text{(Merokok = 1–24)} \\ \lambda^B_3 & \text{(Merokok > 25)} \end{cases} + \lambda^{AB}_{ij} \text{ (jika ada interaksi)} \]
tabel_mi <- matrix(c(25, 25, 12,
0, 1, 3),
nrow = 2, byrow = TRUE)
colnames(tabel_mi) <- c("Merokok0", "Merokok1_24", "Merokok_>25")
rownames(tabel_mi) <- c("Kontrol", "IM")
tabel_mi
## Merokok0 Merokok1_24 Merokok_>25
## Kontrol 25 25 12
## IM 0 1 3
# Ubah menjadi data.frame untuk glm
data_mi <- as.data.frame(as.table(tabel_mi))
colnames(data_mi) <- c("Grup", "Merokok", "Frekuensi")
data_mi
# Model log-linear tanpa interaksi (asumsi independen)
model_tanpa_interaksi <- glm(Frekuensi ~ Grup + Merokok,
family = poisson(),
data = data_mi)
summary(model_tanpa_interaksi)
##
## Call:
## glm(formula = Frekuensi ~ Grup + Merokok, family = poisson(),
## data = data_mi)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.15636 0.20243 15.592 < 2e-16 ***
## GrupIM -2.74084 0.51588 -5.313 1.08e-07 ***
## MerokokMerokok1_24 0.03922 0.28011 0.140 0.889
## MerokokMerokok_>25 -0.51083 0.32660 -1.564 0.118
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 71.5948 on 5 degrees of freedom
## Residual deviance: 6.6901 on 2 degrees of freedom
## AIC: 34.145
##
## Number of Fisher Scoring iterations: 6
# Model log-linear dengan interaksi (untuk cek asosiasi)
model_dengan_interaksi <- glm(Frekuensi ~ Grup * Merokok,
family = poisson(),
data = data_mi)
summary(model_dengan_interaksi)
##
## Call:
## glm(formula = Frekuensi ~ Grup * Merokok, family = poisson(),
## data = data_mi)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.219e+00 2.000e-01 16.094 <2e-16 ***
## GrupIM -2.552e+01 4.225e+04 -0.001 0.9995
## MerokokMerokok1_24 1.402e-15 2.828e-01 0.000 1.0000
## MerokokMerokok_>25 -7.340e-01 3.512e-01 -2.090 0.0366 *
## GrupIM:MerokokMerokok1_24 2.230e+01 4.225e+04 0.001 0.9996
## GrupIM:MerokokMerokok_>25 2.414e+01 4.225e+04 0.001 0.9995
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 7.1595e+01 on 5 degrees of freedom
## Residual deviance: 4.1224e-10 on 0 degrees of freedom
## AIC: 31.455
##
## Number of Fisher Scoring iterations: 20
Karena interaksi tidak signifikan, maka distribusi kebiasaan merokok serupa antara kelompok MI dan kontrol, atau dengan kata lain tidak ada bukti asosiasi kuat antara MI dan merokok dalam kategori ini.
Pada pembahasan sebelumnya, kita telah memahami bahwa salah satu tujuan utama dari penyusunan model log linear adalah untuk mengestimasi parameter-parameter yang menjelaskan hubungan di antara variabel-variabel kategorik.
Pada materi kali ini, kita akan membahas model log linear yang lebih kompleks, yaitu model log linear untuk tabel kontingensi tiga arah. Model ini melibatkan tiga variabel kategorik, sehingga kemungkinan interaksi yang dapat terjadi di dalam model pun menjadi lebih banyak. Dalam konteks ini, interaksi paling tinggi yang dapat dimodelkan adalah interaksi tiga arah, yaitu interaksi yang melibatkan ketiga variabel secara bersamaan.
Model log-linear yang melibatkan tiga variabel kategorik (misal: X, Y, dan Z) dapat dibangun dalam berbagai bentuk model, tergantung pada tingkat interaksi yang ingin dimasukkan. Berikut adalah beberapa alternatif model log-linear yang umum digunakan.
1. Model Saturated
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} + \lambda^{XYZ}_{ijk} \]
Model ini memuat semua kemungkinan interaksi, termasuk interaksi tiga arah (X, Y, dan Z).
2. Model Homogen
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
Model ini hanya mengakomodasi interaksi dua arah antar variabel tanpa memasukkan interaksi tiga arah.
3. Model Conditional
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \]
Memuat interaksi X dengan Y dan X dengan Z.
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \]
Memuat interaksi Y dengan X dan Y dengan Z.
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
Memuat interaksi Z dengan X dan Z dengan Y.
4. Model Joint Independence
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \]
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \]
5. Model Tanpa Interaksi
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k \] Model ini hanya memasukkan efek utama tanpa interaksi antar variabel.
Dalam analisis model log-linear tiga arah, pengujian interaksi dilakukan untuk mengetahui ada atau tidaknya interaksi antar variabel. Pengujian ini dilakukan secara bertahap, dimulai dari tingkat interaksi tertinggi ke yang lebih rendah. Untuk model log-linear dengan tiga peubah (X, Y, dan Z), tahapan pengujian meliputi:
1) Pengujian Interaksi Tiga Arah (XYZ):
2) Pengujian Interaksi Dua Arah (XY, XZ, YZ):
Bandingkan model homogenous dengan model conditional.
Bandingkan model conditional dengan model joint independence.
Bandingkan model joint independence dengan model tanpa interaksi.
Setiap tahapan pengujian dilakukan untuk menilai kecocokan model dan menentukan struktur interaksi mana yang paling sesuai dengan data yang diamati.
Tabel berikut menyajikan data dari survei mengenai penggunaan alkohol, rokok, dan ganja oleh siswa SMA. Susun dan interpretasikan model log-linear paling sederhana (paling parsimonious) untuk data ini. Jelaskan proses yang Anda lakukan dalam menentukan model terbaik serta asosiasi apa saja yang teridentifikasi. Tunjukkan juga bagaimana nilai yang diprediksi dari model menggambarkan asosiasi tersebut.
Alcohol Use | Cigarette Use | Marijuana Use: Yes | Marijuana Use: No |
---|---|---|---|
Yes | Yes | 911 | 538 |
Yes | No | 44 | 456 |
No | Yes | 3 | 43 |
No | No | 2 | 279 |
Keterangan:
library("epitools")
library("DescTools")
library("lawstat")
## Warning: package 'lawstat' was built under R version 4.4.3
Input Data
# Input data sesuai tabel penggunaan Alkohol, Rokok, dan Ganja
alcohol <- factor(rep(c("1Yes", "2No"), each = 4)) # 2 levels × 4 = 8
cigarette <- factor(rep(c("1Yes", "1Yes", "2No", "2No"), 2)) # Ulangin karena nested
marijuana <- factor(rep(c("1Yes", "2No"), times = 4)) # Yes, No → berulang
counts <- c(911, 538, 44, 456, 3, 43, 2, 279)
data <- data.frame(
Alcohol = alcohol,
Cigarette = cigarette,
Marijuana = marijuana,
Frekuensi = counts
)
data
Membentuk Tabel Kontigensi 3 Arah
# Buat table 3 dimensi dari data alkohol, rokok, ganja
table3d <- xtabs(Frekuensi ~ Alcohol + Cigarette + Marijuana, data = data)
# Tampilkan dengan format flat table
ftable(table3d)
## Marijuana 1Yes 2No
## Alcohol Cigarette
## 1Yes 1Yes 911 538
## 2No 44 456
## 2No 1Yes 3 43
## 2No 2 279
Analisis Log-Linear: Tahap Pemodelan Kita akan memodelkan tabel ini menggunakan beberapa model log-linear dan membandingkan kecocokan model (parsimonious model)
data$Alcohol <- relevel(data$Alcohol, ref = "2No")
data$Cigarette <- relevel(data$Cigarette, ref = "2No")
data$Marijuana <- relevel(data$Marijuana, ref = "2No")
Model Saturated Model log-linear saturated memasukkan semua interaksi hingga tiga arah:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} + \lambda^{XYZ}_{ijk} \]
# Model saturated
model_saturated <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana +
Alcohol:Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_saturated)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Alcohol:Marijuana + Cigarette:Marijuana + Alcohol:Cigarette:Marijuana,
## family = poisson(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.63121 0.05987 94.060 < 2e-16
## Alcohol1Yes 0.49128 0.07601 6.464 1.02e-10
## Cigarette1Yes -1.87001 0.16383 -11.414 < 2e-16
## Marijuana1Yes -4.93806 0.70964 -6.959 3.44e-12
## Alcohol1Yes:Cigarette1Yes 2.03538 0.17576 11.580 < 2e-16
## Alcohol1Yes:Marijuana1Yes 2.59976 0.72698 3.576 0.000349
## Cigarette1Yes:Marijuana1Yes 2.27548 0.92746 2.453 0.014149
## Alcohol1Yes:Cigarette1Yes:Marijuana1Yes 0.58951 0.94236 0.626 0.531600
##
## (Intercept) ***
## Alcohol1Yes ***
## Cigarette1Yes ***
## Marijuana1Yes ***
## Alcohol1Yes:Cigarette1Yes ***
## Alcohol1Yes:Marijuana1Yes ***
## Cigarette1Yes:Marijuana1Yes *
## Alcohol1Yes:Cigarette1Yes:Marijuana1Yes
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2.8515e+03 on 7 degrees of freedom
## Residual deviance: 1.3056e-13 on 0 degrees of freedom
## AIC: 65.043
##
## Number of Fisher Scoring iterations: 3
exp(model_saturated$coefficients)
## (Intercept) Alcohol1Yes
## 2.790000e+02 1.634409e+00
## Cigarette1Yes Marijuana1Yes
## 1.541219e-01 7.168459e-03
## Alcohol1Yes:Cigarette1Yes Alcohol1Yes:Marijuana1Yes
## 7.655141e+00 1.346053e+01
## Cigarette1Yes:Marijuana1Yes Alcohol1Yes:Cigarette1Yes:Marijuana1Yes
## 9.732558e+00 1.803106e+00
Model ini adalah model log-linear saturated dengan semua efek utama, interaksi dua arah, dan interaksi tiga arah. Model ini memodelkan hubungan antara Penggunaan Alkohol, Penggunaan Rokok, dan Penggunaan Ganja terhadap frekuensi responden.
library(knitr)
signif_stars <- function(pval) {
if (pval < 0.001) return("***")
else if (pval < 0.01) return("**")
else if (pval < 0.05) return("*")
else if (pval < 0.1) return(".")
else return("")
}
coef_saturated <- round(summary(model_saturated)$coefficients, 3)
colnames(coef_saturated) <- c("Est", "SE", "z", "p")
Signif <- sapply(coef_saturated[, "p"], signif_stars)
coef_table <- cbind(coef_saturated, Sig = Signif)
kable(coef_table, format = "pipe", caption = "Tabel Koefisien Model Saturated")
Est | SE | z | p | Sig | |
---|---|---|---|---|---|
(Intercept) | 5.631 | 0.06 | 94.06 | 0 | *** |
Alcohol1Yes | 0.491 | 0.076 | 6.464 | 0 | *** |
Cigarette1Yes | -1.87 | 0.164 | -11.414 | 0 | *** |
Marijuana1Yes | -4.938 | 0.71 | -6.959 | 0 | *** |
Alcohol1Yes:Cigarette1Yes | 2.035 | 0.176 | 11.58 | 0 | *** |
Alcohol1Yes:Marijuana1Yes | 2.6 | 0.727 | 3.576 | 0 | *** |
Cigarette1Yes:Marijuana1Yes | 2.275 | 0.927 | 2.453 | 0.014 | * |
Alcohol1Yes:Cigarette1Yes:Marijuana1Yes | 0.59 | 0.942 | 0.626 | 0.532 |
(Intercept): Rata-rata log jumlah kasus untuk kategori referensi (tidak minum alkohol, tidak merokok, tidak pakai ganja) adalah 5.631, atau rata-rata kasus sekitar exp(5.631) ≈ 279.
Alcohol1Yes: Mereka yang mengonsumsi alkohol memiliki jumlah kejadian yang 1.63 kali lebih tinggi dibanding yang tidak (exp(0.491) ≈ 1.63), signifikan (p < 0.001).
Cigarette1Yes: Mereka yang merokok memiliki jumlah kejadian sekitar 0.15 kali lebih rendah dari yang tidak merokok (exp(-1.870) ≈ 0.15), signifikan (p < 0.001).
Marijuana1Yes: Mereka yang menggunakan ganja memiliki jumlah kejadian sekitar 0.007 kali lebih rendah dibanding yang tidak (exp(-4.938) ≈ 0.007), signifikan (p < 0.001).
Residual deviance ≈ 0 menandakan model saturated benar-benar fit terhadap data (seluruh variasi data dijelaskan oleh model).
Konsumsi alkohol, konsumsi rokok, dan konsumsi ganja masing-masing secara signifikan meningkatkan jumlah kejadian.
Interaksi antara alkohol dan rokok, antara alkohol dan ganja, serta rokok dan ganja secara signifikan meningkatkan jumlah kejadian.
Interaksi tiga zat (alkohol, rokok, dan ganja) tidak menunjukkan pengaruh yang signifikan.
Model log-linear homogenous memasukkan semua efek utama dan semua interaksi dua arah, tanpa interaksi tiga arah. Secara matematis, model ini dapat dituliskan sebagai berikut:
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]
# Model homogenous (tanpa interaksi tiga arah)
model_homogenous <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_homogenous)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Alcohol:Marijuana + Cigarette:Marijuana, family = poisson(link = "log"),
## data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.63342 0.05970 94.361 < 2e-16 ***
## Alcohol1Yes 0.48772 0.07577 6.437 1.22e-10 ***
## Cigarette1Yes -1.88667 0.16270 -11.596 < 2e-16 ***
## Marijuana1Yes -5.30904 0.47520 -11.172 < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes 2.05453 0.17406 11.803 < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes 2.98601 0.46468 6.426 1.31e-10 ***
## Cigarette1Yes:Marijuana1Yes 2.84789 0.16384 17.382 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.46098 on 7 degrees of freedom
## Residual deviance: 0.37399 on 1 degrees of freedom
## AIC: 63.417
##
## Number of Fisher Scoring iterations: 4
Pengujian ini menggunakan residual deviance dari kedua model (saturated dan homogenous).
H0: Tidak ada interaksi tiga arah (model homogenous sudah cukup)
H1: Ada interaksi tiga arah (model saturated diperlukan)
# Deviance antar model
Deviance.model <- model_homogenous$deviance - model_saturated$deviance
Deviance.model
## [1] 0.3739859
# Derajat bebas = db model homogenous - db model saturated
derajat.bebas <- (model_homogenous$df.residual - model_saturated$df.residual)
derajat.bebas
## [1] 1
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel
## [1] 3.841459
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan
## [1] "Terima H0"
Interpretasi Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara Penggunaan Alkohol, Penggunaan Rokok, dan Penggunaan Ganja. Atau dengan kata lain, model homogenous saja sudah cukup.
Rangkuman
Pengujian Ada Tidaknya Interaksi Tiga Arah (Saturated Model vs Homogenous Model)
Hipotesis:
Tingkat Signifikansi:
Statistik Uji:
\[ \Delta \text{Deviance} = \text{Deviance model homogenous} - \text{Deviance model saturated} \] \[ = 0.374 - 0.0 = 0.374 \]
\[ df = df\ \text{model homogenous} - df\ \text{model saturated} = 1 - 0 = 1 \]
Daerah Penolakan:
Keputusan:
Interpretasi:
Model log-linear conditional pada X memasukkan efek utama dan interaksi dua arah antara X dengan Y dan X dengan Z, tanpa interaksi antara Y dengan Z maupun interaksi tiga arah
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} \]
# Model Conditional Association on X
model_conditional_X <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Alcohol:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_conditional_X)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Alcohol:Marijuana, family = poisson(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.62295 0.06005 93.635 <2e-16 ***
## Alcohol1Yes -0.08167 0.07810 -1.046 0.296
## Cigarette1Yes -1.80971 0.15905 -11.378 <2e-16 ***
## Marijuana1Yes -4.16511 0.45067 -9.242 <2e-16 ***
## Alcohol1Yes:Cigarette1Yes 2.87373 0.16730 17.178 <2e-16 ***
## Alcohol1Yes:Marijuana1Yes 4.12509 0.45294 9.107 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.46 on 7 degrees of freedom
## Residual deviance: 497.37 on 2 degrees of freedom
## AIC: 558.41
##
## Number of Fisher Scoring iterations: 5
• Hipotesis \[ H_0 : \lambda_{jk}^{YZ} = 0 \] Tidak ada interaksi antara Penggunaan Rokok (y) dan Penggunaan Ganja (z)
\[ H_1 : \lambda_{jk}^{YZ} \ne 0 \] Ada interaksi antara Penggunaan Rokok (y) dan Penggunaan Ganja (z)
• Taraf Signifikansi \[ \alpha = 5\% \]
• Statistik Uji \[ \Delta \text{Deviance} = \text{Deviance model conditional on } X - \text{Deviance model homogenous} \]
\[ = 497.369 - 0.374 = 496.995 \]
\[ db = \text{df model conditional on } X - \text{df model homogenous} \]
\[ = 2 - 1 = 1 \]
• Daerah Penolakan \[ \text{Tolak } H_0 \text{ jika } \Delta \text{Deviance} > \chi^2_{0.05, 1} = 3.841 \] • Keputusan
Karena 496.95 > 3.841, maka tolak \(H_0\)
• Kesimpulan Dengan taraf nyata 5%, ada interaksi yang signifikan antara Penggunaan Rokok (Y) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{jk}^{YZ}\).
# Selisih deviance antar model
Deviance.model <- model_conditional_X$deviance - model_homogenous$deviance
Deviance.model
## [1] 496.9953
Hitung Derajat Bebas
# Chi Square tabel dengan alpha = 0.05
derajat.bebas <- (2 - 1)
derajat.bebas
## [1] 1
Nilai Chi-Square Tabel
chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 3.841459
Keputusan Uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Tolak"
Interpretasi
Karena nilai Deviance.model = 996.95 lebih besar dari nilai kritis chi-square tabel = 3.84 (dengan df = 1, alpha = 0.05), maka keputusan uji adalah “Tolak”.
Pada taraf nyata 5%, terdapat bukti yang cukup untuk menolak \(H_0\). Artinya, ada interaksi yang signifkan antara Penggunaan Rokok (Y) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model yang menyertakan parameter interaksi \(\lambda_{jk}^{YZ}\).
Model log-linear conditional pada Y memasukkan efek utama dan interaksi dua arah antara X dengan Y dan Y dengan Z, tanpa interaksi antara X dengan Z maupun interaksi tiga arah.
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{jk}^{YZ} \]
# Model Conditional Association on Y (Cigarette)
model_conditional_Y <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_conditional_Y)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Cigarette:Marijuana, family = poisson(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.57765 0.06032 92.463 < 2e-16 ***
## Alcohol1Yes 0.57625 0.07456 7.729 1.08e-14 ***
## Cigarette1Yes -2.69414 0.16257 -16.572 < 2e-16 ***
## Marijuana1Yes -2.77123 0.15199 -18.233 < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes 2.87373 0.16730 17.178 < 2e-16 ***
## Cigarette1Yes:Marijuana1Yes 3.22431 0.16098 20.029 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.461 on 7 degrees of freedom
## Residual deviance: 92.018 on 2 degrees of freedom
## AIC: 153.06
##
## Number of Fisher Scoring iterations: 6
• Hipotesis
\[ H_0 : \lambda_{ik}^{XZ} = 0 \] Tidak ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Ganja (z)
\[ H_1 : \lambda_{ik}^{XZ} \ne 0 \] Ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Ganja (z)
• Tingkat Signifikansi
\[ \alpha = 5\% \]
• Statistik Uji
\[ \Delta \text{Deviance} = \text{Deviance model conditional on } Y - \text{Deviance model homogenous} \]
\[ = 92.018 - 0.372 = 91.644 \]
\[ db = \text{df model conditional on } Y - \text{df model homogenous} \]
\[ = 2 - 1 = 1 \]
• Daerah Penolakan
Tolak \(H_0\) jika \(\Delta\text{Deviance} > \chi^2_{0.05, 1} = 3.841\)
• Keputusan
\[ \text{Karena } 91.644 > 3.841, \text{ maka tolak } H_0 \] • Kesimpulan
Dengan taraf nyata 5%, ada interaksi yang signifikan antara Penggunaan Alkohol (X) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{ik}^{XZ}\).
# Deviance of Model
Deviance.model <- model_conditional_Y$deviance - model_homogenous$deviance
Deviance.model
## [1] 91.64437
Hitung Derajat Bebas
derajat.bebas <- (2 - 1)
derajat.bebas
## [1] 1
Nilai Chi-Square Tabel
chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 3.841459
Keputusan Uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Tolak"
Interpretasi
Karena nilai Deviance.model = 91.644 lebih besar dari nilai kritis chi-square tabel = 3.84 (dengan df = 1, alpha = 0.05), maka keputusan uji adalah “Tolak”.
Pada taraf nyata 5%, terdapat bukti yang cukup untuk menolak \(H_0\). Artinya, ada interaksi yang signifkan antara Penggunaan Alkohol (X)) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model yang menyertakan parameter interaksi \(\lambda_{ik}^{XZ}\).
Model log-linear conditional pada Z memasukkan efek utama dan interaksi dua arah antara X dengan Z dan Y dengan Z, tanpa interaksi antara X dengan Y maupun interaksi tiga arah. \[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]
# Model Conditional Association on Z (Marijuana)
model_conditional_Z <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Marijuana + Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_conditional_Z)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Marijuana +
## Cigarette:Marijuana, family = poisson(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.19207 0.06088 85.285 < 2e-16 ***
## Alcohol1Yes 1.12719 0.06412 17.579 < 2e-16 ***
## Cigarette1Yes -0.23512 0.05551 -4.235 2.28e-05 ***
## Marijuana1Yes -6.62092 0.47370 -13.977 < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes 4.12509 0.45294 9.107 < 2e-16 ***
## Cigarette1Yes:Marijuana1Yes 3.22431 0.16098 20.029 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.46 on 7 degrees of freedom
## Residual deviance: 187.75 on 2 degrees of freedom
## AIC: 248.8
##
## Number of Fisher Scoring iterations: 5
• Hipotesis
\[ H_0 : \lambda_{ij}^{XY} = 0 \] Tidak ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Rokok (y)
\[ H_1 : \lambda_{ij}^{XY} \ne 0 \] Ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Rokok (y)
• Tingkat Signifikansi
\[ \alpha = 5\% \]
• Statistik Uji
\[ \Delta\text{Deviance} = \text{Deviance model conditional on } Z - \text{Deviance model homogenous} \]
\[ = 187.75 - 0.374 = 187.380 \]
\[ db = \text{df model conditional on } Z - \text{df model homogenous} \]
\[ = 2 - 1 = 1 \]
• Daerah Penolakan
Tolak \(H_0\) jika \(\Delta\text{Deviance} > \chi^2_{0.05, 1} = 3.841\)
• Keputusan
\[ \text{Karena } 187.380 > 3.841, \text{ maka tolak } H_0 \] • Kesimpulan
Dengan taraf nyata 5%, Ada interaksi yang signifikan antara Penggunaan Alkohol (X) dan Penggunaan Rokok (Y)). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{ij}^{XY}\).
# Deviance of Model
Deviance.model <- model_conditional_Z$deviance - model_homogenous$deviance
Deviance.model
## [1] 187.3803
Hitung Derajat Bebas
derajat.bebas <- (2 - 1)
derajat.bebas
## [1] 1
Nilai Chi-Square Tabel
chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 3.841459
Keputusan Uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Tolak"
Interpretasi
Karena nilai Deviance.model = 187.376 lebih besar dari nilai kritis chi-square tabel = 3.84 (dengan df = 1, alpha = 0.05), maka keputusan uji adalah “Tolak”.
Pada taraf nyata 5%, terdapat bukti yang cukup untuk menolak \(H_0\). Artinya, ada interaksi yang signifkan antara Penggunaan Alkohol (X)) dan Penggunaan Rokok (Y)). Dengan kata lain, model yang terbentuk adalah model yang menyertakan parameter interaksi \(\lambda_{ij}^{XY}\).
Model log linier dengan berbagai kombinasi parameter:
library(knitr)
model_table <- data.frame(
Model = c(
"Saturated",
"Homogenous",
"Conditional on X",
"Conditional on Y",
"Conditional on Z"
),
Parameter = c(
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ",
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ",
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ",
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ",
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ"
),
Deviance = c(0.000, 0.374, 497.37, 92.018, 187.75),
Param = c(8, 7, 6, 6, 6),
df = c(0, 1, 1, 1, 1),
AIC = c(65.043, 63.417, 558.41, 153.06, 248.80)
)
kable(
model_table,
format = "pipe",
caption = "Ringkasan Model Log-Linear"
)
Model | Parameter | Deviance | Param | df | AIC |
---|---|---|---|---|---|
Saturated | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ | 0.000 | 8 | 0 | 65.043 |
Homogenous | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ | 0.374 | 7 | 1 | 63.417 |
Conditional on X | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ | 497.370 | 6 | 1 | 558.410 |
Conditional on Y | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ | 92.018 | 6 | 1 | 153.060 |
Conditional on Z | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ | 187.750 | 6 | 1 | 248.800 |
library(knitr)
uji_interaksi <- data.frame(
Interaksi = c("XYZ", "YZ", "XZ", "XY"),
Pengujian = c(
"Saturated vs Homogenous",
"Conditional on X vs Homogenous",
"Conditional on Y vs Homogenous",
"Conditional on Z vs Homogenous"
),
"Δ Deviance" = c(0.374, 496.995, 91.644, 187.380),
df = c(2, 1, 1, 1),
ChiSqTabel = c(5.991, 3.841, 3.841, 3.841),
Ket = c("Tidak ada interaksi", "Ada interaksi", "Ada interaksi", "Ada interaksi"),
Keputusan = c("Terima $H_0$", "Tolak $H_0$", "Tolak $H_0$", "Tolak $H_0$")
)
# Tampilkan dengan kable
kable(uji_interaksi, format = "pipe", caption = "Ringkasan Uji Interaksi (Chi-Square)")
Interaksi | Pengujian | Δ.Deviance | df | ChiSqTabel | Ket | Keputusan |
---|---|---|---|---|---|---|
XYZ | Saturated vs Homogenous | 0.374 | 2 | 5.991 | Tidak ada interaksi | Terima \(H_0\) |
YZ | Conditional on X vs Homogenous | 496.995 | 1 | 3.841 | Ada interaksi | Tolak \(H_0\) |
XZ | Conditional on Y vs Homogenous | 91.644 | 1 | 3.841 | Ada interaksi | Tolak \(H_0\) |
XY | Conditional on Z vs Homogenous | 187.380 | 1 | 3.841 | Ada interaksi | Tolak \(H_0\) |
Dari hasil di atas diketahui bahwa asosiasi yang nyata terdapat pada model Homogenous dan ketiga model Conditional. Karena itulah, model akan dipilih berdasarkan nilai AIC yang paling kecil, dimana AIC terkecil terdapat pada model Homogenous yaitu sebesar 63.417. Oleh karena itu, model terbaik yaitu:
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]
Model terbaik adalah model log-linear tanpa interaksi tiga arah dan memuat semua interaksi dua arah.
Model terbaik dipilih berdasarkan nilai AIC paling kecil, yaitu model Homogenous.
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]
# Model Terbaik
best_model <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(best_model)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Alcohol:Marijuana + Cigarette:Marijuana, family = poisson(link = "log"),
## data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.63342 0.05970 94.361 < 2e-16 ***
## Alcohol1Yes 0.48772 0.07577 6.437 1.22e-10 ***
## Cigarette1Yes -1.88667 0.16270 -11.596 < 2e-16 ***
## Marijuana1Yes -5.30904 0.47520 -11.172 < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes 2.05453 0.17406 11.803 < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes 2.98601 0.46468 6.426 1.31e-10 ***
## Cigarette1Yes:Marijuana1Yes 2.84789 0.16384 17.382 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.46098 on 7 degrees of freedom
## Residual deviance: 0.37399 on 1 degrees of freedom
## AIC: 63.417
##
## Number of Fisher Scoring iterations: 4
Dari summary model diatas terlihat bahwa best model memiliki AIC yang lebih rendah dibandingkan saturated dan conditional model.
# Interpretasi koefisien model terbaik
data.frame(
koef = best_model$coefficients,
exp_koef = exp(best_model$coefficients)
)
\[ \exp(\lambda^{X}_{1A}) = \exp(0.488) = 1.629 \Rightarrow \text{nilai odds} \]
Tanpa memperhatikan penggunaan rokok dan ganja, mereka yang mengonsumsi alkohol memiliki jumlah kasus 1.63 kali lebih besar dibanding yang tidak mengonsumsi alkohol
\[ \exp(\lambda^{Y}_{1C}) = \exp(-1.887) = 0.152 \Rightarrow \text{nilai odds} \]
Tanpa memperhatikan penggunaan alkohol dan ganja, mereka yang merokok memiliki jumlah kasus 0.15 kali lebih kecil dibanding yang tidak merokok
\[ \exp(\lambda^{Z}_{1M}) = \exp(-5.309) = 0.004 \Rightarrow \text{nilai odds} \]
Tanpa memperhatikan penggunaan alkohol dan rokok, mereka yang menggunakan ganja memiliki jumlah kasus hanya sekitar 0.004 kali dari yang tidak menggunakan ganja
\[ \exp(\lambda^{Z}_{1A,1C}) = \exp(2.055) = 7.807 \Rightarrow \text{nilai odds} \]
Responden yang menggunakan Alkohol dan Rokok memiliki jumlah kasus 7.8 kali lebih besar dibanding jika hanya salah satu atau tidak sama sekali.
\[ \exp(\lambda^{XY}_{1A,1M}) = \exp(2.986) = 19.806 \Rightarrow \text{nilai odds ratio} \]
Responden yang menggunakan Alkohol dan Ganja memiliki jumlah kasus 19.8 kali lebih besar dibanding jika hanya salah satu atau tidak sama sekali.
\[ \exp(\lambda^{XY}_{1C, 1M}) = \exp(2.848) = 17.253 \Rightarrow \text{nilai odds ratio} \] Responden yang menggunakan Rokok dan Ganja memiliki jumlah kasus 17.25 kali lebih besar dibanding jika hanya salah satu atau tidak sama sekali.
data.frame(
Alcohol = data$Alcohol,
Cigarette = data$Cigarette,
Marijuana = data$Marijuana,
Frekuensi = data$Frekuensi,
fitted = fitted(model_homogenous)
)
Secara manual, nilai fitted value diperoleh dengan cara sebagai berikut:
\[ \hat{\mu}_{111} = \exp(5.631 + 0.491 - 1.870 - 4.938 + 2.035 + 2.600 + 2.275 + 0.590) = \]
\[ \exp(6.811) = 910.383 \]
\[ \hat{\mu}_{112} = \exp(5.631 + 0.491 - 1.870 + 0 + 2.035 + 0 + 0 + 0) = \]
\[ \exp(6.287) = 538.617 \]
\[ \hat{\mu}_{121} = \exp(5.631 + 0.491 + 0 - 4.938 + 0 + 2.600 + 0 + 0) = \]
\[ \exp(3.799) = 44.617 \]
\[ \hat{\mu}_{122} = \exp(5.631 + 0.491 + 0 + 0 + 0 + 0 + 0 + 0) = \]
\[ \exp(6.122) = 455.383 \]
\[ \hat{\mu}_{211} = \exp(5.631 + 0 + -1.870 - 4.938 + 0 + 0 + 2.275 + 0) = \]
\[ \exp(1.285) = 3.617 \]
\[ \hat{\mu}_{212} = \exp(5.631 + 0 - 1.870 + 0 + 0 + 0 + 0 + 0) = \]
\[ \exp(3.746) = 42.383 \]
\[ \hat{\mu}_{221} = \exp(5.631 + 0 + 0 - 4.938 + 0 + 0 + 0 + 0) = \]
\[ \exp(1.285) = 1.383 \]
\[ \hat{\mu}_{222} = \exp(5.631 + 0 + 0 + 0 + 0 + 0 + 0 + 0) = \]
\[ \exp(5.631) = 279.617 \]
Keterangan: Nilai \(\hat{\mu}_{ijk}\) akan sama apapun referensi dari kategori peubahnya yang kita gunakan
Interpretasi:
Nilai fitted values merupakan prediksi dari model yang sudah mempertimbangkan interaksi dua arah. Nilai ini diperoleh dari estimasi model log-linear yang menyertakan semua efek utama dan semua interaksi dua arah. Dari perbandingan antara prediksi model dan nilai asli, didapatkan bahwa nilai tersebut saling mirip serta tidak terdapat perbedaan yang signifikan. Oleh karena itu, dapat disimpulkan bahwa model dapat memprediksi data dengan sangat baik.
Misalkan kita memiliki tabel tentang pengaruh penggunaan sabuk pengaman terhadap tingkat fatalitas kecelakaan.
Penggunaan Sabuk Pengaman | Cedera Fatal | Cedera Tidak Fatal | Total |
---|---|---|---|
Tidak memakai sabuk pengaman | 1,601 | 162,527 | 164,128 |
Memakai sabuk pengaman | 510 | 412,368 | 412,878 |
Total | 2,111 | 574,895 | 577,006 |
Sumber : Department of Highway Safety and Motor Vehicles, Florida 1988
Hitunglah peluang bersama, peluang marjinal, dan peluang bersyarat beserta uji independensi dan ukuran asosiasi
Peluang Langkah 1 : Hitung Peluang Bersama
Langkah 2 : Hitung Peluang Marginal
Langkah 3 : Hitung Peluang Bersyarat
\[ P(\text{Fatal}|\text{No Seatbelt}) = \frac{1601}{164128} = 0.0098 \]
\[ P(\text{Fatal}|\text{Seatbelt}) = \frac{510}{412878} = 0.0012 \]
\[ P(\text{Seatbelt}|\text{Fatal}) = \frac{510}{2111} = 0.2416 \]
\[ P(\text{Seatbelt}|\text{NonFatal}) = \frac{412368}{574895} = 0.7173 \]
# Membuat tabel 2x2
data <- matrix(c(1601, 162527, 510, 412368), nrow = 2, byrow = TRUE)
rownames(data) <- c("Tidak Memakai Sabuk", "Memakai Sabuk")
colnames(data) <- c("Cedera Fatal", "Cedera Tidak Fatal")
# Total keseluruhan
total <- sum(data)
# Peluang bersama (joint probabilities)
joint_prob <- data / total
joint_prob
## Cedera Fatal Cedera Tidak Fatal
## Tidak Memakai Sabuk 0.002774668 0.2816730
## Memakai Sabuk 0.000883873 0.7146685
# Peluang marginal (baris dan kolom)
marginal_row <- rowSums(data) / total
marginal_col <- colSums(data) / total
marginal_row
## Tidak Memakai Sabuk Memakai Sabuk
## 0.2844476 0.7155524
marginal_col
## Cedera Fatal Cedera Tidak Fatal
## 0.003658541 0.996341459
# Peluang bersyarat
# P(Cedera Fatal | Tidak Memakai Sabuk)
p_fatal_given_no_seatbelt <- data[1, 1] / sum(data[1, ])
# P(Cedera Fatal | Memakai Sabuk)
p_fatal_given_seatbelt <- data[2, 1] / sum(data[2, ])
# P(Memakai Sabuk | Cedera Fatal)
p_seatbelt_given_fatal <- data[2, 1] / sum(data[, 1])
# P(Memakai Sabuk | Cedera Tidak Fatal)
p_seatbelt_given_nonfatal <- data[2, 2] / sum(data[, 2])
p_fatal_given_no_seatbelt
## [1] 0.009754582
p_fatal_given_seatbelt
## [1] 0.001235232
p_seatbelt_given_fatal
## [1] 0.2415917
p_seatbelt_given_nonfatal
## [1] 0.7172927
Interpretasi
Peluang Bersama:
0.00277 → Peluang seseorang tidak memakai sabuk dan cedera fatal.
0.00088 → Peluang seseorang memakai sabuk dan cedera fatal.
0.28167 → Peluang seseorang tidak memakai sabuk dan tidak fatal.
0.71467 → Peluang seseorang memakai sabuk dan tidak fatal.
Peluang Marginal:
0.2844 → Peluang seseorang tidak memakai sabuk secara keseluruhan.
0.7156 → Peluang seseorang memakai sabuk secara keseluruhan.
0.00366 → Peluang mengalami cedera fatal secara umum.
0.9963 → Peluang mengalami cedera tidak fatal secara umum.
Peluang Bersyarat:
Dikarenakan \[ \mathit{P}(\text{Fatal} \mid \text{No Seatbelt}) > \mathit{P}(\text{Fatal} \mid \text{Seatbelt}) \] maka dapat disimpulkan bahwa penggunaan sabuk pengaman dapat menurunkan risiko fatalitas dalam kecelakaan.
Uji Independensi
# Membuat data tabel 2x2
data <- matrix(c(1601, 162527, 510, 412368), nrow = 2, byrow = TRUE)
# Menambahkan nama baris dan kolom
rownames(data) <- c("Tidak Memakai Sabuk", "Memakai Sabuk")
colnames(data) <- c("Cedera Fatal", "Cedera Tidak Fatal")
# Menampilkan data tabel
data
## Cedera Fatal Cedera Tidak Fatal
## Tidak Memakai Sabuk 1601 162527
## Memakai Sabuk 510 412368
# Menghitung uji chi-square
chi_square_test <- chisq.test(data)
# Menampilkan hasil uji chi-square
chi_square_test
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: data
## X-squared = 2336.1, df = 1, p-value < 2.2e-16
Interpretasi:
Dengan taraf signifikansi 5%, didapatkan p-value < 0.05, maka kita menolak H0 (hipotesis nol). Artinya, ada hubungan yang signifikan secara statistik antara penggunaan sabuk pengaman dan tingkat fatalitas kecelakaan
Ukuran Asosiasi
Risk Difference \[ RD = \frac{a}{a + b} - \frac{c}{c + d} = \frac{1601}{164128} - \frac{510}{412878} = 0.0098 - 0.0012 = 0.0086 \]
Relative Risk \[ \frac{P(\text{Fatal} \mid \text{No Seatbelt})}{P(\text{Fatal} \mid \text{Seatbelt})} = \frac{\frac{1601}{164128}}{\frac{510}{412878}} = \frac{0.009759}{0.001235} \approx 7.897 \]
Odds Ratio \[ OR = \frac{a \cdot d}{b \cdot c} = \frac{1601 \cdot 412368}{162527 \cdot 510} = \frac{660201168}{82988770} = 7.96 \]
# Data
a <- 1601
b <- 162527
c <- 510
d <- 412368
# Risk Difference
RD <- (a / (a + b)) - (c / (c + d))
# Relative Risk
RR <- (a / (a + b)) / (c / (c + d))
# Odds Ratio
OR <- (a * d) / (b * c)
#Hasil
list(Risk_Difference = RD, Relative_Risk = RR, Odds_Ratio = OR)
## $Risk_Difference
## [1] 0.00851935
##
## $Relative_Risk
## [1] 7.896965
##
## $Odds_Ratio
## [1] 7.964905
Interpretasi :
RD = 0.0086 → Perbedaan risiko kecelakaan fatal sebesar 0.86% lebih tinggi pada pengendara yang tidak pakai sabuk.
RR = 7.90 → Pengendara yang tidak memakai sabuk memiliki 8.17 kali risiko lebih besar dalam mengalami kecelakaan fatal.
OR = 7.96 → Peluang kecelakaan fatal 7.96 kali lebih besar pada pengendara yang tidak memakai sabuk.
Tabel berikut menyajikan data dari survei mengenai penggunaan alkohol, rokok, dan ganja oleh mahasiswa. Tentukan model log-linear terbaik yang akan digunakan dan tunjukkan bagaimana nilai yang diprediksi dari model menggambarkan asosiasi pada data tersebut
Alcohol Use | Cigarette Use | Marijuana Use: Yes | Marijuana Use: No |
---|---|---|---|
Yes | Yes | 911 | 538 |
Yes | No | 44 | 456 |
No | Yes | 3 | 43 |
No | No | 2 | 279 |
Sumber: Wright State University School of Medicine & United Health Services. (1992). Survey of substance use among high school seniors in a nonurban area near Dayton, Ohio. Dayton, OH
Keterangan:
Data
library("epitools")
library("DescTools")
library("lawstat")
# Input data sesuai tabel penggunaan Alkohol, Rokok, dan Ganja
alcohol <- factor(rep(c("1Yes", "2No"), each = 4)) # 2 levels × 4 = 8
cigarette <- factor(rep(c("1Yes", "1Yes", "2No", "2No"), 2)) # Ulangin karena nested
marijuana <- factor(rep(c("1Yes", "2No"), times = 4)) # Yes, No → berulang
counts <- c(911, 538, 44, 456, 3, 43, 2, 279)
data <- data.frame(
Alcohol = alcohol,
Cigarette = cigarette,
Marijuana = marijuana,
Frekuensi = counts
)
#Relevel Data
data$Alcohol <- relevel(data$Alcohol, ref = "2No")
data$Cigarette <- relevel(data$Cigarette, ref = "2No")
data$Marijuana <- relevel(data$Marijuana, ref = "2No")
# Buat table 3 dimensi dari data alkohol, rokok, ganja
table3d <- xtabs(Frekuensi ~ Alcohol + Cigarette + Marijuana, data = data)
# Tampilkan dengan format flat table
ftable(table3d)
## Marijuana 2No 1Yes
## Alcohol Cigarette
## 2No 2No 279 2
## 1Yes 43 3
## 1Yes 2No 456 44
## 1Yes 538 911
Perbandingan Model
Model Saturated
Model log-linear saturated memasukkan semua interaksi hingga tiga arah:
\[ \log(\mu_{ijk}) = \lambda + \lambda^X_i + \lambda^Y_j + \lambda^Z_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} + \lambda^{XYZ}_{ijk} \]
# Model saturated
model_saturated <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana +
Alcohol:Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_saturated)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Alcohol:Marijuana + Cigarette:Marijuana + Alcohol:Cigarette:Marijuana,
## family = poisson(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.63121 0.05987 94.060 < 2e-16
## Alcohol1Yes 0.49128 0.07601 6.464 1.02e-10
## Cigarette1Yes -1.87001 0.16383 -11.414 < 2e-16
## Marijuana1Yes -4.93806 0.70964 -6.959 3.44e-12
## Alcohol1Yes:Cigarette1Yes 2.03538 0.17576 11.580 < 2e-16
## Alcohol1Yes:Marijuana1Yes 2.59976 0.72698 3.576 0.000349
## Cigarette1Yes:Marijuana1Yes 2.27548 0.92746 2.453 0.014149
## Alcohol1Yes:Cigarette1Yes:Marijuana1Yes 0.58951 0.94236 0.626 0.531600
##
## (Intercept) ***
## Alcohol1Yes ***
## Cigarette1Yes ***
## Marijuana1Yes ***
## Alcohol1Yes:Cigarette1Yes ***
## Alcohol1Yes:Marijuana1Yes ***
## Cigarette1Yes:Marijuana1Yes *
## Alcohol1Yes:Cigarette1Yes:Marijuana1Yes
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2.8515e+03 on 7 degrees of freedom
## Residual deviance: 1.3056e-13 on 0 degrees of freedom
## AIC: 65.043
##
## Number of Fisher Scoring iterations: 3
Model Homogenous
Model log-linear homogenous memasukkan semua efek utama dan semua interaksi dua arah, tanpa interaksi tiga arah. Secara matematis, model ini dapat dituliskan sebagai berikut:
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]
# Model homogenous
model_homogenous <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Alcohol:Marijuana + Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_homogenous)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Alcohol:Marijuana + Cigarette:Marijuana, family = poisson(link = "log"),
## data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.63342 0.05970 94.361 < 2e-16 ***
## Alcohol1Yes 0.48772 0.07577 6.437 1.22e-10 ***
## Cigarette1Yes -1.88667 0.16270 -11.596 < 2e-16 ***
## Marijuana1Yes -5.30904 0.47520 -11.172 < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes 2.05453 0.17406 11.803 < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes 2.98601 0.46468 6.426 1.31e-10 ***
## Cigarette1Yes:Marijuana1Yes 2.84789 0.16384 17.382 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.46098 on 7 degrees of freedom
## Residual deviance: 0.37399 on 1 degrees of freedom
## AIC: 63.417
##
## Number of Fisher Scoring iterations: 4
Uji Hipotesis: Apakah Ada Interaksi Tiga Arah (Saturated vs Homogenous)
Pengujian ini menggunakan residual deviance dari kedua model (saturated dan homogenous).
Hipotesis
H0: Tidak ada interaksi tiga arah (model homogenous sudah cukup)
H1: Ada interaksi tiga arah (model saturated diperlukan)
Selisih Deviance
# Deviance antar model
Deviance.model <- model_homogenous$deviance - model_saturated$deviance
Deviance.model
## [1] 0.3739859
Derajat Bebas
# Derajat bebas = db model homogenous - db model saturated
derajat.bebas <- (model_homogenous$df.residual - model_saturated$df.residual)
derajat.bebas
## [1] 1
Chi Square Tabel
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel
## [1] 3.841459
Keputusan Uji
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan
## [1] "Terima H0"
Interpretasi
Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara Penggunaan Alkohol, Penggunaan Rokok, dan Penggunaan Ganja. Atau dengan kata lain, model homogenous saja sudah cukup.
Model Conditional on X
Model log-linear conditional pada X memasukkan efek utama dan interaksi dua arah antara X dengan Y dan X dengan Z, tanpa interaksi antara Y dengan Z maupun interaksi tiga arah
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} \]
# Model Conditional Association on X
model_conditional_X <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Alcohol:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_conditional_X)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Alcohol:Marijuana, family = poisson(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.62295 0.06005 93.635 <2e-16 ***
## Alcohol1Yes -0.08167 0.07810 -1.046 0.296
## Cigarette1Yes -1.80971 0.15905 -11.378 <2e-16 ***
## Marijuana1Yes -4.16511 0.45067 -9.242 <2e-16 ***
## Alcohol1Yes:Cigarette1Yes 2.87373 0.16730 17.178 <2e-16 ***
## Alcohol1Yes:Marijuana1Yes 4.12509 0.45294 9.107 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.46 on 7 degrees of freedom
## Residual deviance: 497.37 on 2 degrees of freedom
## AIC: 558.41
##
## Number of Fisher Scoring iterations: 5
Pengujian Ada Tidaknya Interaksi Antara Y dan Z (Homogenous Model vs Conditional Association on X)
• Hipotesis \[ H_0 : \lambda_{jk}^{YZ} = 0 \] Tidak ada interaksi antara Penggunaan Rokok (y) dan Penggunaan Ganja (z)
\[ H_1 : \lambda_{jk}^{YZ} \ne 0 \] Ada interaksi antara Penggunaan Rokok (y) dan Penggunaan Ganja (z)
• Taraf Signifikansi \[ \alpha = 5\% \]
• Statistik Uji \[ \Delta \text{Deviance} = \text{Deviance model conditional on } X - \text{Deviance model homogenous} \]
\[ = 497.369 - 0.374 = 496.995 \]
\[ db = \text{df model conditional on } X - \text{df model homogenous} \]
\[ = 2 - 1 = 1 \]
• Daerah Penolakan \[ \text{Tolak } H_0 \text{ jika } \Delta \text{Deviance} > \chi^2_{0.05, 1} = 3.841 \] • Keputusan
Karena 496.95 > 3.841, maka tolak \(H_0\)
• Kesimpulan Dengan taraf nyata 5%, ada interaksi yang signifikan antara Penggunaan Rokok (Y) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{jk}^{YZ}\).
Model Conditional on Y
Model log-linear conditional pada Y memasukkan efek utama dan interaksi dua arah antara X dengan Y dan Y dengan Z, tanpa interaksi antara X dengan Z maupun interaksi tiga arah.
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{jk}^{YZ} \]
# Model Conditional Association on Y (Cigarette)
model_conditional_Y <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Cigarette + Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_conditional_Y)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Cigarette +
## Cigarette:Marijuana, family = poisson(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.57765 0.06032 92.463 < 2e-16 ***
## Alcohol1Yes 0.57625 0.07456 7.729 1.08e-14 ***
## Cigarette1Yes -2.69414 0.16257 -16.572 < 2e-16 ***
## Marijuana1Yes -2.77123 0.15199 -18.233 < 2e-16 ***
## Alcohol1Yes:Cigarette1Yes 2.87373 0.16730 17.178 < 2e-16 ***
## Cigarette1Yes:Marijuana1Yes 3.22431 0.16098 20.029 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.461 on 7 degrees of freedom
## Residual deviance: 92.018 on 2 degrees of freedom
## AIC: 153.06
##
## Number of Fisher Scoring iterations: 6
Pengujian Ada Tidaknya Interaksi antara X dan Z (Homogenous model vs Conditional Association on Y)
• Hipotesis
\[ H_0 : \lambda_{ik}^{XZ} = 0 \] Tidak ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Ganja (z)
\[ H_1 : \lambda_{ik}^{XZ} \ne 0 \] Ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Ganja (z)
• Tingkat Signifikansi
\[ \alpha = 5\% \]
• Statistik Uji
\[ \Delta \text{Deviance} = \text{Deviance model conditional on } Y - \text{Deviance model homogenous} \]
\[ = 92.018 - 0.372 = 91.644 \]
\[ db = \text{df model conditional on } Y - \text{df model homogenous} \]
\[ = 2 - 1 = 1 \]
• Daerah Penolakan
Tolak \(H_0\) jika \(\Delta\text{Deviance} > \chi^2_{0.05, 1} = 3.841\)
• Keputusan
\[ \text{Karena } 91.644 > 3.841, \text{ maka tolak } H_0 \] • Kesimpulan
Dengan taraf nyata 5%, ada interaksi yang signifikan antara Penggunaan Alkohol (X) dan Penggunaan Ganja (Z). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{ik}^{XZ}\).
Model Conditional on Z
Model log-linear conditional pada Z memasukkan efek utama dan interaksi dua arah antara X dengan Z dan Y dengan Z, tanpa interaksi antara X dengan Y maupun interaksi tiga arah. \[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]
# Model Conditional Association on Z (Marijuana)
model_conditional_Z <- glm(counts ~ Alcohol + Cigarette + Marijuana +
Alcohol:Marijuana + Cigarette:Marijuana,
family = poisson(link = "log"), data = data)
summary(model_conditional_Z)
##
## Call:
## glm(formula = counts ~ Alcohol + Cigarette + Marijuana + Alcohol:Marijuana +
## Cigarette:Marijuana, family = poisson(link = "log"), data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 5.19207 0.06088 85.285 < 2e-16 ***
## Alcohol1Yes 1.12719 0.06412 17.579 < 2e-16 ***
## Cigarette1Yes -0.23512 0.05551 -4.235 2.28e-05 ***
## Marijuana1Yes -6.62092 0.47370 -13.977 < 2e-16 ***
## Alcohol1Yes:Marijuana1Yes 4.12509 0.45294 9.107 < 2e-16 ***
## Cigarette1Yes:Marijuana1Yes 3.22431 0.16098 20.029 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2851.46 on 7 degrees of freedom
## Residual deviance: 187.75 on 2 degrees of freedom
## AIC: 248.8
##
## Number of Fisher Scoring iterations: 5
Pengujian Ada Tidaknya Interaksi antara X dan Y (Homogenous model vs Conditional Association on Z)
• Hipotesis
\[ H_0 : \lambda_{ij}^{XY} = 0 \] Tidak ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Rokok (y)
\[ H_1 : \lambda_{ij}^{XY} \ne 0 \] Ada interaksi antara Penggunaan Alkohol (x) dan Penggunaan Rokok (y)
• Tingkat Signifikansi
\[ \alpha = 5\% \]
• Statistik Uji
\[ \Delta\text{Deviance} = \text{Deviance model conditional on } Z - \text{Deviance model homogenous} \]
\[ = 187.75 - 0.374 = 187.380 \]
\[ db = \text{df model conditional on } Z - \text{df model homogenous} \]
\[ = 2 - 1 = 1 \]
• Daerah Penolakan
Tolak \(H_0\) jika \(\Delta\text{Deviance} > \chi^2_{0.05, 1} = 3.841\)
• Keputusan
\[ \text{Karena } 187.380 > 3.841, \text{ maka tolak } H_0 \] • Kesimpulan
Dengan taraf nyata 5%, Ada interaksi yang signifikan antara Penggunaan Alkohol (X) dan Penggunaan Rokok (Y)). Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{ij}^{XY}\).
Pemilihan Model Terbaik
Ringkasan Model Log Linier
Model log linier dengan berbagai kombinasi parameter:
library(knitr)
model_table <- data.frame(
Model = c(
"Saturated",
"Homogenous",
"Conditional on X",
"Conditional on Y",
"Conditional on Z"
),
Parameter = c(
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ",
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ",
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ",
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ",
"λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ"
),
Deviance = c(0.000, 0.374, 497.37, 92.018, 187.75),
Param = c(8, 7, 6, 6, 6),
df = c(0, 1, 1, 1, 1),
AIC = c(65.043, 63.417, 558.41, 153.06, 248.80)
)
kable(
model_table,
format = "pipe",
caption = "Ringkasan Model Log-Linear"
)
Model | Parameter | Deviance | Param | df | AIC |
---|---|---|---|---|---|
Saturated | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ | 0.000 | 8 | 0 | 65.043 |
Homogenous | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ | 0.374 | 7 | 1 | 63.417 |
Conditional on X | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ | 497.370 | 6 | 1 | 558.410 |
Conditional on Y | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ | 92.018 | 6 | 1 | 153.060 |
Conditional on Z | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ | 187.750 | 6 | 1 | 248.800 |
Ringkasan Pengujian Interaksi 3 Arah dan 2 Arah
library(knitr)
uji_interaksi <- data.frame(
Interaksi = c("XYZ", "YZ", "XZ", "XY"),
Pengujian = c(
"Saturated vs Homogenous",
"Conditional on X vs Homogenous",
"Conditional on Y vs Homogenous",
"Conditional on Z vs Homogenous"
),
"Δ Deviance" = c(0.374, 496.995, 91.644, 187.380),
df = c(2, 1, 1, 1),
ChiSqTabel = c(5.991, 3.841, 3.841, 3.841),
Ket = c("Tidak ada interaksi", "Ada interaksi", "Ada interaksi", "Ada interaksi"),
Keputusan = c("Terima $H_0$", "Tolak $H_0$", "Tolak $H_0$", "Tolak $H_0$")
)
# Tampilkan dengan kable
kable(uji_interaksi, format = "pipe", caption = "Ringkasan Uji Interaksi (Chi-Square)")
Interaksi | Pengujian | Δ.Deviance | df | ChiSqTabel | Ket | Keputusan |
---|---|---|---|---|---|---|
XYZ | Saturated vs Homogenous | 0.374 | 2 | 5.991 | Tidak ada interaksi | Terima \(H_0\) |
YZ | Conditional on X vs Homogenous | 496.995 | 1 | 3.841 | Ada interaksi | Tolak \(H_0\) |
XZ | Conditional on Y vs Homogenous | 91.644 | 1 | 3.841 | Ada interaksi | Tolak \(H_0\) |
XY | Conditional on Z vs Homogenous | 187.380 | 1 | 3.841 | Ada interaksi | Tolak \(H_0\) |
Kesimpulan Pemilihan Model Terbaik
Dari hasil di atas diketahui bahwa asosiasi yang nyata terdapat pada model Homogenous dan ketiga model Conditional. Karena itulah, model akan dipilih berdasarkan nilai AIC yang paling kecil, dimana AIC terkecil terdapat pada model Homogenous yaitu sebesar 63.417. Oleh karena itu, model terbaik yaitu:
\[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]
Model terbaik adalah model log-linear tanpa interaksi tiga arah dan memuat semua interaksi dua arah.
Nilai Dugaan Model Terbaik
data.frame(
Alcohol = data$Alcohol,
Cigarette = data$Cigarette,
Marijuana = data$Marijuana,
Frekuensi = data$Frekuensi,
fitted = fitted(model_homogenous)
)
Interpretasi:
Nilai fitted values merupakan prediksi dari model yang sudah mempertimbangkan interaksi dua arah. Nilai ini diperoleh dari estimasi model log-linear yang menyertakan semua efek utama dan semua interaksi dua arah. Dari perbandingan antara prediksi model dan nilai asli, didapatkan bahwa nilai tersebut saling mirip serta tidak terdapat perbedaan yang signifikan. Oleh karena itu, dapat disimpulkan bahwa model dapat memprediksi data dengan sangat baik.
Agresti, A. (2002). Categorical data analysis (2nd ed.). Wiley-Interscience.
Agresti, A. (2013). Categorical data analysis (3rd ed.). John Wiley & Sons.
Agresti, A. (2019). An introduction to categorical data analysis (3rd ed.). Wiley.
Christensen, R. (1997). Log-linear models and logistic regression. Springer.
Fisher, R. A. (1935). The design of experiments. Oliver and Boyd.
Fox, J. (2008). Applied regression analysis and generalized linear models (2nd ed.). Sage.
Hilbe, J. M. (2009). Logistic regression models. Chapman and Hall/CRC.
Hosmer, D. W., Lemeshow, S., & Sturdivant, R. X. (2013). Applied logistic regression (3rd ed.). Wiley.
Howell, D. C. (2012). Statistical methods for psychology (8th ed.). Cengage Learning.
Kleinbaum, D. G., & Klein, M. (2010). Logistic regression: A self-learning text (3rd ed.). Springer.
Mehta, C. R., & Patel, N. R. (1983). A network algorithm for performing Fisher’s exact test in r × c contingency tables. Journal of the American Statistical Association, 78(382), 427–434.
Modul ADK, Pertemuan 13. (2025). Universitas Padjadjaran.
Powers, D. M. W. (2011). Evaluation: From precision, recall and F-measure to ROC, informedness, markedness and correlation. Journal of Machine Learning Technologies, 2(1), 37–63