Dalam analisis data, tabel kontingensi digunakan untuk memahami hubungan antara dua variabel kategorikal. Tabel kontingensi adalah tabel yang menyajikan distribusi frekuensi dari dua variabel dalam bentuk matriks. Tabel ini digunakan untuk menganalisis hubungan antara variabel dan menguji signifikansi statistik hubungan tersebut.
Tabel kontingensi 2x2 adalah bentuk paling sederhana dari tabel kontingensi, di mana terdapat dua kategori untuk masing-masing variabel. Dengan tabel ini, kita dapat menghitung berbagai ukuran probabilitas dan asosiasi antara variabel yang diamati.
Pada contoh ini, kita akan menganalisis hubungan antara kebiasaan olahraga dan kondisi kesehatan menggunakan tabel kontingensi 2x2.
Bentuk Umum Tabel Kontingensi 2x2:
Tabel kontingensi 2x2 biasanya disusun dalam format berikut:
\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Kategori 1} & \textbf{Kategori 2} & \textbf{Total} \\ \hline \textbf{Kelompok 1} & n_{11} & n_{12} & n_{1.} \\ \textbf{Kelompok 2} & n_{21} & n_{22} & n_{2.} \\ \hline \textbf{Total} & n_{.1} & n_{.2} & n \\ \hline \end{array} \]
Dimana:
\(n_{11}\): Jumlah kasus kategori 1 dari kelompok 1
\(n_{12}\): Jumlah kasus kategori 2 dari kelompok 1
\(n_{21}\): Jumlah kasus kategori 2 dari kelompok 1
\(n_{22}\): Jumlah kasus kategori 2 dari kelompok 2
\(n_{1.}\): Total kasus dari kelompok 1
\(n_{2.}\): Total kasus dari kelompok 2
\(n_{.1}\): Total kasus dari kategori 1
\(n_{.2}\): Total kasus dari kategori 2
\(n\): Total jumlah kasus keseluruhan
Dalam tabel kontingensi 2x2 terdapat beberapa distribusi peluang, yaitu peluang bersama, peluang marginal, dan peluang bersyarat.
Peluang bersama adalah probabilitas bahwa dua kejadian terjadi secara bersamaan. Dalam tabel kontingensi 2x2, peluang bersama dihitung sebagai:
\[ P(A_i \cap B_j) = \frac{n_{ij}}{n} \]
Dimana:
\(n_{ij}\) adalah jumlah observasi pada sel (i, j).
\(n\) adalah total jumlah observasi.
Peluang marginal adalah probabilitas kejadian tunggal tanpa memperhatikan kejadian lain. Peluang marginal untuk suatu baris atau kolom dihitung sebagai:
\[ P(A_i) = \frac{n_{i+}}{n}, \quad P(B_j) = \frac{n_{+j}}{n} \]
Dimana:
\(n_{i+}\) adalah total observasi pada baris i.
\(n_{+j}\) adalah total observasi pada kolom j.
\(n\) adalah total jumlah observasi.
Peluang bersyarat adalah probabilitas suatu kejadian terjadi dengan syarat kejadian lain sudah terjadi. Dalam tabel kontingensi, peluang bersyarat dihitung sebagai:
\[ P(A_i | B_j) = \frac{P(A_i \cap B_j)}{P(B_j)} = \frac{n_{ij} / n}{n_{+j} / n} = \frac{n_{ij}}{n_{+j}} \]
atau
\[ P(B_j | A_i) = \frac{n_{ij}}{n_{i+}} \]
Sebuah tim peneliti dari sebuah universitas kesehatan masyarakat melakukan survei kepada 100 orang dewasa di kota X. Survei ini bertujuan untuk mengetahui apakah terdapat hubungan antara kebiasaan berolahraga dengan kondisi kesehatan seseorang. Dalam survei ini, responden diminta untuk mengisi kuesioner yang mencakup frekuensi mereka dalam melakukan aktivitas olahraga dan penilaian umum terhadap kondisi kesehatan mereka.
Peneliti mengelompokkan kebiasaan olahraga responden menjadi dua kategori: - Berolahraga - Tidak Berolahraga
Sementara itu, kondisi kesehatan dikelompokkan berdasarkan penilaian subyektif responden terhadap kesehatannya menjadi dua kategori: - Sehat - Tidak Sehat
Hasil dari survei ini disajikan dalam bentuk tabel kontingensi berikut: \[ \begin{array}{|c|c|c|c|} \hline & \textbf{Sehat (B1)} & \textbf{Tidak Sehat (B2)} & \textbf{Total} \\ \hline \textbf{Olahraga (A1)} & 30 & 10 & 40 \\ \textbf{Tidak Olahraga (A2)} & 20 & 40 & 60 \\ \hline \textbf{Total} & 50 & 50 & 100 \\ \hline \end{array} \]
Hitung Peluang Bersama
\(P(A_1 \cap B_1) = \frac{30}{100} = 0.3\)
\(P(A_1 \cap B_2) = \frac{10}{100} = 0.1\)
\(P(A_2 \cap B_1) = \frac{20}{100} = 0.2\)
\(P(A_2 \cap B_2) = \frac{40}{100} = 0.4\)
Dengan menggunakan syntax R:
n_ij <- matrix(c(30, 10, 20, 40), nrow = 2, byrow = TRUE)
colnames(n_ij) <- c("Sehat", "Tidak Sehat")
rownames(n_ij) <- c("Olahraga", "Tidak Olahraga")
n <- sum(n_ij)
p_join <- n_ij / n
print(p_join)
## Sehat Tidak Sehat
## Olahraga 0.3 0.1
## Tidak Olahraga 0.2 0.4
Hitung Peluang Marginal
Peluang marginal untuk variabel olahraga:
\(P(A_1) = P(A_1 \cap B_1) + P(A_1 \cap B_2) = 0.3 + 0.1 = 0.4\)
\(P(A_2) = P(A_2 \cap B_1) + P(A_2 \cap B_2) = 0.2 + 0.4 = 0.6\)
Peluang marginal untuk variabel kesehatan:
n_ij <- matrix(c(30, 10, 20, 40), nrow = 2, byrow = TRUE)
p_marginal_A <- rowSums(p_join)
p_marginal_B <- colSums(p_join)
list(Peluang_Marginal_A = p_marginal_A, Peluang_Marginal_B = p_marginal_B)
## $Peluang_Marginal_A
## Olahraga Tidak Olahraga
## 0.4 0.6
##
## $Peluang_Marginal_B
## Sehat Tidak Sehat
## 0.5 0.5
Hitung Peluang Bersyarat
\[ P(B1 | A1) = \frac{P(A1 \cap B1)}{P(A1)} = \frac{30}{40} = 0.75 \]
\[ P(B2 | A1) = \frac{P(A1 \cap B2)}{P(A1)} = \frac{10}{40} = 0.25 \]
\[ P(B1 | A2) = \frac{P(A2 \cap B1)}{P(A2)} = \frac{20}{60} = 0.33 \]
\[ P(B2 | A2) = \frac{P(A2 \cap B2)}{P(A2)} = \frac{40}{60} = 0.67 \]
n_ij <- matrix(c(30, 10, 20, 40), nrow = 2, byrow = TRUE)
colnames(n_ij) <- c("Sehat", "Tidak Sehat")
rownames(n_ij) <- c("Olahraga", "Tidak Olahraga")
p_conditional <- n_ij / rowSums(n_ij)
print(p_conditional)
## Sehat Tidak Sehat
## Olahraga 0.7500000 0.2500000
## Tidak Olahraga 0.3333333 0.6666667Interpretasi:
Peluang bersama menunjukkan peluang gabungan dari kejadian tertentu dalam tabel kontingensi 2x2
Peluang marginal menunjukkan peluang terjadinya suatu kejadian tanpa memperhatikan variabel lain
Peluang bersyarat menunjukkan peluang terjadinya suatu kejadian berubah ketika memperhatikan variabel lain
Dikarenakan, \(P(Sehat|Olahraga)=0.75 > P(Sehat|Tidak Olahraga)=0.25\), maka menunjukkan bahwa berolahraga lebih memberikan dampak sehat dibandingkan tidak berolahraga.
Ukuran asosiasi dalam tabel kontingensi 2x2 digunakan untuk mengukur hubungan antara dua variabel. Ukuran asosiasi yang umum digunakan, yaitu Risk Difference (RD)/Beda Peluang, Relative Risk/Resiko Relatif (RR), dan Odds Ratio (OR).
Risk Difference (RD) atau beda peluang adalah selisih antara probabilitas kejadian dalam dua kelompok. RD menunjukkan seberapa besar perbedaan peluang antara kelompok yang terpapar dan kelompok yang tidak terpapar.
\[ RD = P(B1 | A1) - P(B1 | A2) \]
Jika RD = 0, berarti tidak terdapat perbedaan peluang antara dua kelompok tersebut.
Relative Risk (RR) atau risiko relatif adalah perbandingan probabilitas kejadian dalam dua kelompok. RR digunakan untuk mengetahui seberapa besar risiko kejadian dalam kelompok terpapar dibandingkan dengan kelompok tidak terpapar.
\[ RR = \frac{P(B1 | A1)}{P(B1 | A2)} \]
Jika RR = 1, berarti tidak terdapat perbedaan risiko antara kedua kelompok tersebut.
Odds Ratio (OR) adalah ukuran yang membandingkan odds kejadian suatu peristiwa dalam dua kelompok. OR sering digunakan dalam studi kasus-kontrol.
\[ OR = \frac{(n_{11}/n_{12})}{(n_{21}/n_{22})} = \frac{n_{11} \cdot n_{22}}{n_{12} \cdot n_{21}} \] Dimana:
\(n_{11}\) = jumlah individu yang mengalami kejadian dalam kelompok pertama
\(n_{12}\) = jumlah individu yang tidak mengalami kejadian dalam kelompok pertama
\(n_{21}\) = jumlah individu yang mengalami kejadian dalam kelompok kedua
\(n_{22}\) = jumlah individu yang tidak mengalami kejadian dalam kelompok kedua
Jika OR = 1, tidak terdapat asosiasi antara 2 kelompok tersebut.
Masih dengan kasus yang sama, kali ini akan dilihat lebih jauh seperti hubungan antara olahraga dengan kesehatan seseorang.
Hitung Risk Difference (RD)
\[ RD = \frac{30}{40} - \frac{20}{60} = 0.75 - 0.33 = 0.42 \]
RD <- (n_ij[1,1] / sum(n_ij[1,])) - (n_ij[2,1] / sum(n_ij[2,]))
print(RD)
## [1] 0.4166667Hitung Relative Risk (RR)
\[ RR = \frac{30/40}{20/60} = \frac{0.75}{0.33} = 2.27 \]
RR <- (n_ij[1,1] / sum(n_ij[1,])) / (n_ij[2,1] / sum(n_ij[2,]))
print(RR)
## [1] 2.25Hitung Odds Ratio (OR)
\[ OR = \frac{(30 \times 40)}{(10 \times 20)} = \frac{1200}{200} = 6 \]
OR <- (n_ij[1,1] * n_ij[2,2]) / (n_ij[1,2] * n_ij[2,1])
print(OR)
## [1] 6Interpretasi:
RD = 0.42, maka orang yang berolahraga memiliki peluang 0.42 lebih tinggi untuk sehat dibandingkan dengan orang yang tidak berolahraga.
RR = 2.27, maka orang yang berolahraga 2.27 kali lebih mungkin untuk sehat dibandingkan dengan orang yang tidak berolahraga.
OR = 6, maka orang yang berolahraga 6 kali lebih mungkin untuk sehat dibandingkan dengan orang yang tidak berolahraga.
Dalam statistik, inferensi adalah suatu proses pengampilan kesimpulan mengenai populasi berdasarkan sampel data. Dalam konteks tabel kontingensi dua arah, inferensi dilakukan untuk menganalisis hubungan antara dua variabel kategorik yang ada dalam tabel kontingensi. Terdapat 2 jenis inferensi dalam tabel kontingensi, yaitu estimasi dan pengujian.
Estimasi statistik adalah proses untuk memperkirakan nilai suatu parameter populasi berdasarkan data dari sampel. Karena data populasi sering tidak tersedia secara keseluruhan, kita menggunakan data sampel untuk memperkirakan parameter populasi seperti rata-rata (mean), proporsi (proportion), atau variansi (variance).
Estimasi terbagi menjadi dua jenis utama:
Estimasi titik memberikan perkiraan tunggal terhadap suatu parameter populasi.
Contoh:
Notasi:
\[ \bar{x} = \frac{1}{n} \sum_{i=1}^{n} x_i \]
\[ \hat{p} = \frac{x}{n} \]
di mana: - \(\bar{x}\): rata-rata sampel - \(\hat{p}\): proporsi sampel - \(x_i\): data ke-i - \(n\): ukuran sampel - \(x\): jumlah kejadian “sukses”
Estimasi interval memberikan rentang nilai di sekitar estimasi titik, yang diduga kuat mencakup nilai parameter populasi dengan tingkat kepercayaan tertentu, seperti 95%.
Jika data berdistribusi normal dan simpangan baku populasi tidak diketahui, maka digunakan distribusi t-Student:
\[ \bar{x} \pm t_{\alpha/2, \, df} \cdot \frac{s}{\sqrt{n}} \]
Keterangan: - \(\bar{x}\): rata-rata
sampel
- \(s\): simpangan baku sampel
- \(n\): ukuran sampel
- \(df = n - 1\): derajat bebas
- \(t_{\alpha/2}\): nilai t pada
tingkat kepercayaan tertentu
\[ \hat{p} \pm z_{\alpha/2} \cdot \sqrt{ \frac{ \hat{p}(1 - \hat{p}) }{n} } \]
Keterangan: - \(\hat{p}\): proporsi sampel - \(z_{\alpha/2}\): nilai z dari distribusi normal standar - \(n\): ukuran sampel
Sebuah kelompok peneliti di bidang kesehatan masyarakat ingin mengetahui apakah konsumsi vitamin C secara rutin dapat mengurangi risiko terkena flu. Untuk itu, mereka melakukan studi terhadap 200 orang dewasa selama musim flu berlangsung (sekitar 3 bulan).
Responden dibagi menjadi dua kelompok berdasarkan kebiasaan mereka mengonsumsi suplemen vitamin C:
Setelah musim flu berakhir, para peneliti mencatat apakah masing-masing responden terkena flu atau tidak terkena flu selama periode pengamatan. Data hasil pengamatan disusun dalam tabel kontingensi 2x2 berikut:
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Terkena Flu (B1)} & \textbf{Tidak Terkena Flu (B2)} & \textbf{Total} \\ \hline \textbf{Minum Vitamin (A1)} & 25 & 75 & 100 \\ \textbf{Tidak Minum Vitamin (A2)} & 50 & 50 & 100 \\ \hline \textbf{Total} & 75 & 125 & 200 \\ \hline \end{array} \]
Uji proporsi dalam tabel kontingensi 2x2 digunakan untuk membandingkan proporsi dua kelompok dan menentukan apakah perbedaan antara proporsi tersebut signifikan secara statistik.
Hipotesis:
Statistik Uji:
\[ Z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{p(1-p)\left(\frac{1}{n_1} + \frac{1}{n_2}\right)}} \] Dimana:
\(\hat{p}_1\) dan \(\hat{p}_2\) adalah proporsi sampel dari masing-masing kelompok.
\(p\) adalah proporsi gabungan dari kedua kelompok.
\(n_1\) dan \(n_2\) adalah ukuran sampel masing-masing kelompok.
Dengan,
\(\hat{p}_1 = \frac{n_{11}}{n_{1.}}\)
\(\hat{p}_2 = \frac{n_{21}}{n_{2.}}\)
\(\hat{p} = \frac{n_{11}+n_{21}}{n_{1.}+n_{2.}}\)
Kriteria Uji:
Jika nilai Z lebih besar dari nilai kritis dalam distribusi normal standar, maka H₀ ditolak, yang berarti terdapat perbedaan signifikan antara dua proporsi.
Jika nilai p-value lebih kecil dari tingkat signifikansi (misalnya 0.05), maka H₀ ditolak.
Hipotesis
Hitung Proporsi Sampel
\[ \hat{p}_1 = \frac{25}{100} = 0.25, \quad \hat{p}_2 = \frac{50}{100} = 0.50 \]
Hitung Proporsi Gabungan
\[ \hat{p} = \frac{25 + 50}{100 + 100} = \frac{75}{200} = 0.375 \]
Hitung Standar Error
\[ SE = \sqrt{p(1-p) \left(\frac{1}{n_1} + \frac{1}{n_2}\right)} \] \[ SE = \sqrt{0.375 \times 0.625 \left(\frac{1}{100} + \frac{1}{100}\right)} \] \[ SE = \sqrt{0.2344 \times 0.02} = \sqrt{0.004688} = 0.0685 \]
Hitung Statistik Uji
\[ Z = \frac{\hat{p}_1 - \hat{p}_2}{SE} \] \[ Z = \frac{0.25 - 0.50}{0.0685} = \frac{-0.25}{0.0685} = -3.65 \]
Tentukan Kriteria Uji
Untuk tingkat signifikansi \(\alpha = 0.05\), nilai kritis dari distribusi normal standar adalah 1.96.
Karena nilai Z = -3.65 lebih kecil dari -1.96, maka H₀ ditolak, yang berarti terdapat perbedaan signifikan antara kelompok yang mengonsumsi vitamin dan yang tidak dalam kemungkinan terkena flu.
Dengan menggunakan syntax R:
data_vitamin <- matrix(c(25, 75, 50, 50), nrow = 2, byrow = TRUE)
rownames(data_vitamin) <- c("Minum Vitamin", "Tidak Minum Vitamin")
colnames(data_vitamin) <- c("Terkena Flu", "Tidak Terkena Flu")
# Uji proporsi untuk data konsumsi vitamin
test_result_vitamin <- prop.test(x = c(25, 50), n = c(100, 100), alternative = "two.sided", correct = FALSE)
# Menampilkan hasil uji proporsi
test_result_vitamin
##
## 2-sample test for equality of proportions without continuity correction
##
## data: c(25, 50) out of c(100, 100)
## X-squared = 13.333, df = 1, p-value = 0.0002607
## alternative hypothesis: two.sided
## 95 percent confidence interval:
## -0.3796394 -0.1203606
## sample estimates:
## prop 1 prop 2
## 0.25 0.50
Interpretasi:
Dari hasil perhitungan di atas, didapatkan bahwa terdapat perbedaan signifikan antara kelompok yang mengonsumsi vitamin C dan yang tidak mengonsumsi vitamin C terhadap kemungkinan terkena flu.
Maka, konsumsi vitamin C memiliki pengaruh signifikan terhadap kemungkinan terkena flu. Seseorang yang mengonsumsi vtamin C memiliki kemungkinan lebih kecil terkena flu dibandingkan dengan yang tidak mengonsumsi vitamin C
Untuk menguji signifikansi dari ketiga ukuran asosiasi, kita dapat menggunakan statistik Z. Asumsinya adalah data berasal dari sampel besar sehingga pendekatan normal dapat digunakan. Dengan hipotesis sebagai berikut:
Hipotesis nol (H₀): Tidak ada asosiasi antara dua kelompok
Hipotesis alternatif (H₁): Ada asosiasi antara dua kelompok
Statistik Uji:
Risk Difference (RD)
Standard error Risk Difference (RD)
\[ SE_{RD} = \sqrt{ \frac{\hat{p}_1(1 - \hat{p}_1)}{n_{1.}} + \frac{\hat{p}_2(1 - \hat{p}_2)}{n_{2.}} } \]
Statistik uji Risk Difference (RD)
\[ Z_{RD} = \frac{RD}{SE_{RD}} \]
Relative Risk (RR)
Standard error Relative Risk (RR)
\[ SE_{\log RR} = \sqrt{ \frac{1}{n_{11}} - \frac{1}{n_{1.}} + \frac{1}{n_{21}} - \frac{1}{n_{2.}} } \]
Statistik Uji Relative Risk (RR)
\[ Z_{RR} = \frac{\log(RR)}{SE_{\log RR}} \]
Odds Ratio (OR)
Standard error Odds Ratio (OR)
\[ SE_{\log OR} = \sqrt{ \frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}} } \]
Statistik Uji Z Odds Ratio
\[ Z_{OR} = \frac{\log(OR)}{SE_{\log OR}} \]
Kriteria uji:
Jika nilai Z lebih besar dari nilai kritis dalam distribusi normal standar, maka H₀ ditolak, yang berarti terdapat asosiasi anatara 2 kelompok.
Jika nilai p-value lebih kecil dari tingkat signifikansi (misalnya 0.05), maka H₀ ditolak.
Contoh kasus:
Misalkan dilakukan penelitian untuk mengetahui apakah konsumsi vitamin C dapat mengurangi kemungkinan seseorang terkena flu. Berikut adalah hasil penelitian dalam tabel kontingensi 2x2:
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Terkena Flu (B1)} & \textbf{Tidak Terkena Flu (B2)} & \textbf{Total} \\ \hline \textbf{Minum Vitamin (A1)} & 25 & 75 & 100 \\ \textbf{Tidak Minum Vitamin (A2)} & 50 & 50 & 100 \\ \hline \textbf{Total} & 75 & 125 & 200 \\ \hline \end{array} \]
Risk Difference (RD)
\[ RD = \hat{p}_1 - \hat{p}_2 = 0.25 - 0.50 = -0.25 \]
\[ SE_{RD} = \sqrt{ \frac{\hat{p}_1(1 - \hat{p}_1)}{n_{1.}} + \frac{\hat{p}_2(1 - \hat{p}_2)}{n_{2.}} } = \sqrt{ \frac{0.25 \times 0.75}{100} + \frac{0.50 \times 0.50}{100} } = \sqrt{ 0.001875 + 0.0025 } = \sqrt{0.004375} \approx 0.0661 \]
\[ Z_{RD} = \frac{RD}{SE_{RD}} = \frac{-0.25}{0.0661} \approx -3.78 \]
Relative Risk (RR)
\[ RR = \frac{\hat{p}_1}{\hat{p}_2} = \frac{0.25}{0.50} = 0.5 \]
\[ SE_{\log RR} = \sqrt{ \frac{1}{n_{11}} - \frac{1}{n_{1.}} + \frac{1}{n_{21}} - \frac{1}{n_{2.}} } = \sqrt{ \frac{1}{25} - \frac{1}{100} + \frac{1}{50} - \frac{1}{100} } = \sqrt{0.04 - 0.01 + 0.02 - 0.01} = \sqrt{0.04} = 0.2 \]
\[ Z_{RR} = \frac{\log(RR)}{SE_{\log RR}} = \frac{\log(0.5)}{0.2} = \frac{-0.6931}{0.2} \approx -3.47 \]
Odds Ratio (OR)
\[ OR = \frac{n_{11} \cdot n_{22}}{n_{12} \cdot n_{21}} = \frac{25 \cdot 50}{75 \cdot 50} = \frac{1250}{3750} = 0.333 \]
\[ SE_{\log OR} = \sqrt{ \frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}} } = \sqrt{ \frac{1}{25} + \frac{1}{75} + \frac{1}{50} + \frac{1}{50} } = \sqrt{0.04 + 0.0133 + 0.02 + 0.02} = \sqrt{0.0933} \approx 0.3055 \]
\[ Z_{OR} = \frac{\log(OR)}{SE_{\log OR}} = \frac{\log(0.333)}{0.3055} = \frac{-1.0986}{0.3055} \approx -3.60 \]
Dengan menggunakan syntax R:
# Data
a <- 25 # Minum vitamin, terkena flu
b <- 75 # Minum vitamin, tidak flu
c <- 50 # Tidak minum, terkena flu
d <- 50 # Tidak minum, tidak flu
n1 <- a + b # Total yang minum
n2 <- c + d # Total yang tidak minum
# Hitung proporsi risiko
p1 <- a / n1 # Risiko kelompok minum
p2 <- c / n2 # Risiko kelompok tidak minum
## 1. Risk Difference (RD)
RD <- p1 - p2
SE_RD <- sqrt( (p1 * (1 - p1)) / n1 + (p2 * (1 - p2)) / n2 )
Z_RD <- RD / SE_RD
p_RD <- 2 * pnorm(-abs(Z_RD))
## 2. Risk Ratio (RR)
RR <- p1 / p2
SE_logRR <- sqrt((1 - p1) / (a) + (1 - p2) / (c))
Z_logRR <- log(RR) / SE_logRR
p_RR <- 2 * pnorm(-abs(Z_logRR))
## 3. Odds Ratio (OR)
OR <- (a * d) / (b * c)
SE_logOR <- sqrt(1/a + 1/b + 1/c + 1/d)
Z_logOR <- log(OR) / SE_logOR
p_OR <- 2 * pnorm(-abs(Z_logOR))
## Tampilkan hasil
cat("=== RISK DIFFERENCE (RD) ===\n")
## === RISK DIFFERENCE (RD) ===
cat("RD:", RD, "\nZ:", Z_RD, "\np-value:", p_RD, "\n\n")
## RD: -0.25
## Z: -3.779645
## p-value: 0.0001570523
cat("=== RISK RATIO (RR) ===\n")
## === RISK RATIO (RR) ===
cat("RR:", RR, "\nZ:", Z_logRR, "\np-value:", p_RR, "\n\n")
## RR: 0.5
## Z: -3.465736
## p-value: 0.0005287824
cat("=== ODDS RATIO (OR) ===\n")
## === ODDS RATIO (OR) ===
cat("OR:", OR, "\nZ:", Z_logOR, "\np-value:", p_OR, "\n")
## OR: 0.3333333
## Z: -3.596053
## p-value: 0.0003230822
Interpretasi:
Untuk tingkat signifikansi \(\alpha = 0.05\), nilai kritis dari distribusi normal standar adalah 1.96.
Karena nilai Z = -3.78, -3.47, -3.60 lebih kecil dari -1.96, maka H₀ ditolak, yang berarti terdapat asosiasi antara 2 kelompok.
Maka, antara kelompok konsumsi vitamin C memiliki asosiasi kelompok kemungkinan terkena flu.
Uji independensi pada tabel kontingensi 2x2 adalah metode statistik untuk menentukan apakah dua variabel kategori saling bebas (independen) atau berhubungan (dependen) satu sama lain. Tujuannya untuk menentukan apakah ada asosiasi atau ketergantungan antara variabel baris (A) dan kolom (B).
Hipotesis
Hipotesis Nol (H₀): Kedua variabel independen (tidak saling berhubungan).
Hipotesis Alternatif (H₁): Kedua variabel tidak independen (ada hubungan/asosiasi).
Statistik Uji:
\[ \chi^2 = \sum \frac{(O - E)^2}{E} \]
Dengan:
Frekuensi harapan dihitung sebagai:
\[ E_{ij} = \frac{(\text{total baris}_i) \cdot (\text{total kolom}_j)}{n} \]
Kriteria Uji:
Jika \(\chi^2_{\text{hitung}} > \chi^2_{\text{tabel}}\), maka tolak H₀
Jika \(\chi^2_{\text{hitung}} \leq \chi^2_{\text{tabel}}\), maka gagal tolak H₀
Jika p-value < 0.05 → Tolak \(H_0\), artinya ada hubungan antara dua variabel.
Jika p-value ≥ 0.05 → Gagal tolak \(H_0\), artinya tidak ada bukti cukup bahwa kedua variabel saling berhubungan.
Contoh Kasus:
Misalkan dilakukan penelitian untuk mengetahui apakah konsumsi vitamin C dapat mengurangi kemungkinan seseorang terkena flu. Berikut adalah hasil penelitian dalam tabel kontingensi 2x2:
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Terkena Flu (B1)} & \textbf{Tidak Terkena Flu (B2)} & \textbf{Total} \\ \hline \textbf{Minum Vitamin (A1)} & 25 & 75 & 100 \\ \textbf{Tidak Minum Vitamin (A2)} & 50 & 50 & 100 \\ \hline \textbf{Total} & 75 & 125 & 200 \\ \hline \end{array} \]
Hipotesis:
H₀ (Hipotesis Nol): Tidak ada hubungan antara minum vitamin dan kejadian flu (variabel independen).
H₁ (Hipotesis Alternatif): Ada hubungan antara minum vitamin dan kejadian flu (variabel tidak independen).
Statistik uji:
Hitung Nilai Frekuensi Harapan:
\[ E_{ij} = \frac{(\text{Total Baris}_i) \times (\text{Total Kolom}_j)}{n} \]
Hitung Statistik uji Chi-Square
\[ \chi^2 = \sum \frac{(O - E)^2}{E} \]
\[\begin{align*} \chi^2 &= \frac{(25 - 37.5)^2}{37.5} + \frac{(75 - 62.5)^2}{62.5} + \frac{(50 - 37.5)^2}{37.5} + \frac{(50 - 62.5)^2}{62.5} \\ &= \frac{(-12.5)^2}{37.5} + \frac{(12.5)^2}{62.5} + \frac{(12.5)^2}{37.5} + \frac{(-12.5)^2}{62.5} \\ &= \frac{156.25}{37.5} + \frac{156.25}{62.5} + \frac{156.25}{37.5} + \frac{156.25}{62.5} \\ &= 4.1667 + 2.5 + 4.1667 + 2.5 \\ &= 13.3334 \end{align*}\]
Hitung nilai Chi-Square Tabel
Derajat kebebasan: \(df = (2 - 1)(2 - 1) =
1\) - Taraf signifikansi: \(\alpha =
0.05\) - Nilai kritis dari tabel Chi-Square:
\[
\chi^2_{(0.05;1)} = 3.841
\]
Dengan menggunakan syntax R:
# Buat tabel kontingensi
flu_matrix <- matrix(c(25, 75, 50, 50), nrow = 2, byrow = TRUE)
# Tambahkan nama baris dan kolom
colnames(flu_matrix) <- c("Flu", "Tidak_Flu")
rownames(flu_matrix) <- c("Minum_Vitamin", "Tidak_Minum")
flu_table <- as.table(flu_matrix)
# Uji chi-square independensi
chi_result <- chisq.test(flu_table, correct = FALSE) # continuity correction dinonaktifkan untuk hasil lebih murni
# Tampilkan hasil
print(chi_result)
##
## Pearson's Chi-squared test
##
## data: flu_table
## X-squared = 13.333, df = 1, p-value = 0.0002607
Keputusan
Untuk tingkat signifikansi \(\alpha = 0.05\), nilai kritis dari distribusi chi-square adalah 3.841.
Karena nilai X2 = 13.3334 lebih besar dari 3.841, maka H₀ ditolak, yang berarti terdapat hubungan antara variabel konsumsi vitamin C dengan kemungkinan terkena flu.
Kesimpulan
Maka, antara kelompok konsumsi vitamin C memiliki hubungan kelompok kemungkinan terkena flu dan tidak saling bebas (tidak independen).
Partisi chi-square adalah prosedur penguraian (dekomposisi) nilai chi-square total dari tabel kontingensi berukuran I × J (baris × kolom) menjadi komponen-komponen chi-square yang lebih kecil. Ini berguna untuk mengetahui sumber utama penyimpangan antara data yang diamati dan yang diharapkan.
Misalnya: jika hasil uji chi-square menunjukkan adanya hubungan antara dua variabel, partisi chi-square dapat membantu menjelaskan dimana tepatnya hubungan tersebut terjadi, yaitu di baris atau kolom tertentu.
Langkah-langkah:
Hitung nilai Chi-Square total
Jika signifikan, maka lanjutkan ke partisi
Tentukan strategi partisi
Hitung nilai chi-square dengan data yang sudah di partisi
Bandingkan dengan nilai kritis.
Contoh Kasus:
Sebuah tim peneliti sosial dari sebuah universitas ingin mengetahui apakah terdapat hubungan antara jenis kelamin dan preferensi minuman favorit masyarakat. Penelitian ini dilakukan dengan menyebarkan kuesioner kepada 150 responden dari berbagai kalangan usia dan profesi.
Dalam kuesioner, responden diminta untuk memilih salah satu jenis minuman yang paling mereka sukai dari tiga pilihan: - Kopi - Teh - Jus
Selain itu, identitas jenis kelamin juga dicatat, dan dikelompokkan ke dalam dua kategori: - Pria - Wanita
Setelah data dikumpulkan dan diklasifikasikan, hasilnya ditampilkan dalam tabel kontingensi berikut:
\[\begin{array}{|c|c|c|c|c|} \hline \textbf{Jenis Kelamin} & \textbf{Kopi} & \textbf{Teh} & \textbf{Jus} & \textbf{Total} \\ \hline \text{Pria} & 30 & 20 & 10 & 60 \\ \text{Wanita} & 10 & 30 & 20 & 60 \\ \hline \text{Total} & 40 & 50 & 30 & 120 \\ \hline \end{array}\]
# Matriks Observasi
observed <- matrix(c(30, 20, 10, 10, 30, 20),
nrow = 2, byrow = TRUE,
dimnames = list(
Jenis_Kelamin = c("Pria", "Wanita"),
Minuman = c("Kopi", "Teh", "Jus")
))
observed
## Minuman
## Jenis_Kelamin Kopi Teh Jus
## Pria 30 20 10
## Wanita 10 30 20
Hitung nilai chi-square total
\[ \begin{array}{|c|c|c|c|} \hline E_{ij} & \text{Kopi} & \text{Teh} & \text{Jus} \\ \hline \text{Pria} & \frac{60 \times 40}{120} = 20 & \frac{60 \times 50}{120} = 25 & \frac{60 \times 30}{120} = 15 \\ \text{Wanita} & 20 & 25 & 15 \\ \hline \end{array} \]
\[\begin{align*} \chi^2 &= \frac{(30 - 20)^2}{20} + \frac{(20 - 25)^2}{25} + \frac{(10 - 15)^2}{15} \\ &\quad + \frac{(10 - 20)^2}{20} + \frac{(30 - 25)^2}{25} + \frac{(20 - 15)^2}{15} \\ &= \frac{100}{20} + \frac{25}{25} + \frac{25}{15} + \frac{100}{20} + \frac{25}{25} + \frac{25}{15} \\ &= 5 + 1 + 1.67 + 5 + 1 + 1.67 \\ &= 15.34 \end{align*}\]
Dengan derajat bebas:
\[ df = (2 - 1)(3 - 1) = 2 \]
Bandingkan nilai \(\chi^2 = 15.34\) dengan nilai kritis \(\chi^2_{0.05, 2} = 5.991\)
Karena \(15.34 > 5.991\), maka terdapat hubungan signifikan antara jenis kelamin dan preferensi minuman.
# Hitung frekuensi harapan
expected <- chisq.test(observed)$expected
expected
## Minuman
## Jenis_Kelamin Kopi Teh Jus
## Pria 20 25 15
## Wanita 20 25 15
# Uji chi-square dari R
chi_test <- chisq.test(observed, correct = FALSE)
chi_test
##
## Pearson's Chi-squared test
##
## data: observed
## X-squared = 15.333, df = 2, p-value = 0.0004682
Melanjutkan ke partisi chi-square
Buat tabel 2x2 baru:
| Kopi+Teh | Jus | Total | |
|---|---|---|---|
| Pria | 50 | 10 | 60 |
| Wanita | 40 | 20 | 60 |
| Total | 90 | 30 | 120 |
# Gabungkan kolom Kopi dan Teh menjadi satu
observed_partisi <- matrix(c(50, 10, 40, 20),
nrow = 2, byrow = TRUE,
dimnames = list(
Jenis_Kelamin = c("Pria", "Wanita"),
Minuman = c("Kopi+Teh", "Jus")
))
observed_partisi
## Minuman
## Jenis_Kelamin Kopi+Teh Jus
## Pria 50 10
## Wanita 40 20
Hitung nilai chi-square partisi:
\[ E_{11} = \frac{60 \times 90}{120} = 45,\quad E_{12} = \frac{60 \times 30}{120} = 15 \]
\[ \chi^2_{\text{partisi}} = \frac{(50 - 45)^2}{45} + \frac{(10 - 15)^2}{15} + \frac{(40 - 45)^2}{45} + \frac{(20 - 15)^2}{15} \]
\[ = \frac{25}{45} + \frac{25}{15} + \frac{25}{45} + \frac{25}{15} = 0.56 + 1.67 + 0.56 + 1.67 = 4.46 \]
expected_partisi <- chisq.test(observed_partisi)$expected
expected_partisi
## Minuman
## Jenis_Kelamin Kopi+Teh Jus
## Pria 45 15
## Wanita 45 15
chi_partisi_test <- chisq.test(observed_partisi, correct = FALSE)
chi_partisi_test
##
## Pearson's Chi-squared test
##
## data: observed_partisi
## X-squared = 4.4444, df = 1, p-value = 0.03501
Dengan \(df = 1\), nilai kritis \(\chi^2_{0.05,1} = 3.841\)
Karena \(4.46 > 3.841\), maka partisi ini signifikan, artinya penyimpangan utama terjadi antara kelompok yang memilih Jus dan lainnya.
Interpretasi:
Uji chi-square menunjukkan hubungan signifikan.
Partisi chi-square mengindikasikan bahwa penyimpangan utama terjadi pada pilihan minuman Jus antara pria dan wanita.
Uji Likelihood Ratio (G2) adalah alternatif dari uji Chi-Square Pearson yang digunakan untuk menguji apakah dua variabel kategori bersifat independen. Dalam konteks tabel kontingensi 2x2, uji ini berguna untuk memverifikasi apakah ada asosiasi antara dua variabel.
Hipotesis
Rumus Statistik Likelihood Ratio
Rumus untuk menghitung nilai statistik Likelihood Ratio (G²) pada tabel kontingensi adalah sebagai berikut:
\[ G^2 = 2 \sum_{i=1}^{r} \sum_{j=1}^{c} O_{ij} \cdot \ln \left( \frac{O_{ij}}{E_{ij}} \right) \]
Di mana:
\[ E_{ij} = \frac{(\text{Total Baris}_i) \cdot (\text{Total Kolom}_j)}{N} \]
Derajat Bebas
Untuk tabel kontingensi berukuran 2x2:
\[ df = (r - 1)(c - 1) = (2 - 1)(2 - 1) = 1 \]
Kriteria Uji:
Jika \(G^2_{\text{hitung}} > \chi^2_{\text{tabel}}\), maka tolak H₀
Jika \(G^2_{\text{hitung}} \leq \chi^2_{\text{tabel}}\), maka gagal tolak H₀
Jika p-value < 0.05 → Tolak \(H_0\), artinya ada hubungan antara dua variabel.
Jika p-value ≥ 0.05 → Gagal tolak \(H_0\), artinya tidak ada bukti cukup bahwa kedua variabel saling berhubungan.
Contoh Kasus:
Sebuah penelitian ingin mengetahui apakah terdapat hubungan antara jenis kelamin dan status merokok.
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Merokok (B1)} & \textbf{Tidak Merokok (B2)} & \textbf{Total} \\ \hline \textbf{Pria (A1)} & 30 & 20 & 50 \\ \textbf{Wanita (A2)} & 10 & 40 & 50 \\ \hline \textbf{Total} & 40 & 60 & 100 \\ \hline \end{array} \]
# Data Observasi
observed <- matrix(c(30, 20, 10, 40),
nrow = 2, byrow = TRUE,
dimnames = list(
Jenis_Kelamin = c("Pria", "Wanita"),
Merokok = c("Merokok", "Tidak_Merokok")
))
observed
## Merokok
## Jenis_Kelamin Merokok Tidak_Merokok
## Pria 30 20
## Wanita 10 40
# Marginal dan total
row_totals <- rowSums(observed)
col_totals <- colSums(observed)
grand_total <- sum(observed)
row_totals
## Pria Wanita
## 50 50
col_totals
## Merokok Tidak_Merokok
## 40 60
grand_total
## [1] 100
Hitung frekuensi ekspektasi
\(E_{11} = \frac{50 \cdot 40}{100} = 20\)
\(E_{12} = \frac{50 \cdot 60}{100} = 30\)
\(E_{21} = \frac{50 \cdot 40}{100} = 20\)
\(E_{22} = \frac{50 \cdot 60}{100} = 30\)
# Hitung nilai E_ij manual
E11 <- (row_totals[1] * col_totals[1]) / grand_total
E12 <- (row_totals[1] * col_totals[2]) / grand_total
E21 <- (row_totals[2] * col_totals[1]) / grand_total
E22 <- (row_totals[2] * col_totals[2]) / grand_total
expected_manual <- matrix(c(E11, E12, E21, E22),
nrow = 2, byrow = TRUE,
dimnames = list(
Jenis_Kelamin = c("Pria", "Wanita"),
Merokok = c("Merokok", "Tidak_Merokok")
))
expected_manual
## Merokok
## Jenis_Kelamin Merokok Tidak_Merokok
## Pria 20 30
## Wanita 20 30
Hitung statistik uji (G2)
Hitung satu per satu:
Jumlah total:
\[ G^2 = 2 \cdot (12.165 - 8.11 - 6.931 + 11.508) = 2 \cdot 8.632 = 17.264 \]
# Hitung G² secara manual
G2_manual <- 2 * (
observed[1,1] * log(observed[1,1] / expected_manual[1,1]) +
observed[1,2] * log(observed[1,2] / expected_manual[1,2]) +
observed[2,1] * log(observed[2,1] / expected_manual[2,1]) +
observed[2,2] * log(observed[2,2] / expected_manual[2,2])
)
G2_manual
## [1] 17.26092
Hitung derajat bebas
\[ df = (r - 1)(c - 1) = (2 - 1)(2 - 1) = 1 \]
p_value_manual <- pchisq(G2_manual, df = 1, lower.tail = FALSE)
p_value_manual
## [1] 3.258188e-05Keputusan
Nilai \(G^2 = 17.264\), dibandingkan dengan nilai kritis \(\chi^2_{0.05,1} = 3.841\). Karena:
\[ G^2 = 17.264 > 3.841 \]
Dan \[ \text{p-value} \approx 0.00003 \]
Karena p-value < 0.05 atau \(G^2\) > nilai kritis, maka kita menolak \(H_0\).
Kesimpulan
Maka, antara kelompok jenis kelamin memiliki hubungan dengan kelompok merokok (tidak independen).
Uji Independensi Eksak Fisher (Fisher’s Exact Test) digunakan untuk menguji apakah dua variabel kategorik pada tabel kontingensi berukuran kecil (umumnya 2x2) saling bebas (independen) atau tidak. Uji ini sangat tepat digunakan ketika frekuensi harapan kecil, khususnya jika ada sel dengan nilai harapan kurang dari 5.
Uji ini tidak mengandalkan distribusi chi-square dan memberikan nilai p eksak, berbeda dari pendekatan aproksimasi.
Hipotesis:
H₀ (Hipotesis nol): Kedua variabel bersifat independen (tidak ada hubungan).
H₁ (Hipotesis alternatif): Kedua variabel tidak independen (ada hubungan).
Rumus Distribusi Hipergeometrik:
Untuk menghitung peluang konfigurasi tabel kontingensi tertentu, digunakan:
\[ P(X = n_{11}) = \frac{{\binom{n_{1.}}{n_{11}} \cdot \binom{n_{2.}}{n_{.1} - n_{11}}}}{{\binom{n}{n_{.1}}}} \]
Keterangan:
Contoh Kasus:
Sebuah studi ingin melihat hubungan antara jenis kelamin dan preferensi minuman:
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Kopi (B1)} & \textbf{Teh (B2)} & \textbf{Total} \\ \hline \textbf{Pria (A1)} & 8 & 2 & 10 \\ \textbf{Wanita (A2)} & 1 & 5 & 6 \\ \hline \textbf{Total} & 40 & 60 & 100 \\ \hline \end{array} \]
Kita ingin menghitung peluang terjadinya tabel yang sama ekstrem atau lebih ekstrem daripada tabel observasi:
\[ P(X = 8) = \frac{\binom{10}{8} \cdot \binom{6}{1}}{\binom{16}{9}} = \frac{45 \cdot 6}{11440} = \frac{270}{11440} \approx 0.0236 \]
Peluang untuk \(X = 9\) (lebih ekstrem mendukung H₁):
\[ P(X = 9) = \frac{\binom{10}{9} \cdot \binom{6}{0}}{\binom{16}{9}} = \frac{10 \cdot 1}{11440} = \frac{10}{11440} \approx 0.00087 \]
Peluang untuk \(X = 3\) (lebih ekstrem mendukung H₁):
\[ P(X = 3) = \frac{\binom{10}{3} \cdot \binom{6}{6}}{\binom{16}{9}} = \frac{120 \cdot 1}{11440} = \frac{10}{11440} \approx 0.01049 \]
Total p-value (dua sisi konservatif):
\[ p\text{-value} = P(X=8) + P(X=9) + P(X=3) = 0.0236 + 0.00087 + 0.01049 = \boxed{0.03496} \]
Jika, menggunakan fisher.test() dengan R:
# Data tabel kontingensi
data <- matrix(c(8, 2, 1, 5),
nrow = 2,
byrow = TRUE,
dimnames = list("Jenis Kelamin" = c("Pria", "Wanita"),
"Minuman" = c("Kopi", "Teh")))
data
## Minuman
## Jenis Kelamin Kopi Teh
## Pria 8 2
## Wanita 1 5
# Uji Fisher
fisher.test(data, alternative = "two.sided")
##
## Fisher's Exact Test for Count Data
##
## data: data
## p-value = 0.03497
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 1.008849 1049.791446
## sample estimates:
## odds ratio
## 15.46969
Kesimpulan:
Dengan p-value sebesar 0.03496 < 0.05, kita menolak H₀. Terdapat bukti bahwa terdapat hubungan antara jenis kelamin dan preferensi minuman dalam sampel tersebut.
Tabel kontingensi 3 arah adalah ekstensi dari tabel kontingensi dua arah yang digunakan untuk menyajikan dan menganalisis hubungan antara tiga variabel kategorik. Dalam bentuk ini, data disusun dalam array tiga dimensi dengan ukuran \(I \times J \times K\), di mana:
Tujuan Analisis:
Analisis tabel kontingensi 3 arah bertujuan untuk:
Struktur Tabel Kontingensi 3 Arah:
Dengan tabel marjinal sebagai berikut:
Total frekuensi pada masing-masing variabel bisa dihitung dengan:
Data menampilkan tabel kontingensi 3 arah, yang memperlihatkan hubungan pengaruh jenis kelamin, konsumsi minuman beralkohol terhadap perilaku seks bebas pada remaja. Data berjumlah 620 dengan proporsi 350 laki-laki, dan 270 perempuan.
\[\begin{array}{|c|c|c|c|} \hline \textbf{Jenis Kelamin} & \textbf{Alkohol} & \textbf{Seks Bebas: Ya} & \textbf{Seks Bebas: Tidak} \\ \hline Laki-Laki & Ya & 130 & 100 \\ Laki-Laki & Tidak & 30 & 90 \\ Perempuan & Ya & 70 & 90 \\ Perempuan & Tidak & 20 & 90 \\ \hline \end{array}\]
Tabel marginal: jumlah total untuk kombinasi dua variabel (X dan Y) dengan mengabaikan yang ketiga (Z). Tabel parsial: tabel dua arah pada kondisi tetap satu variabel (Z) tertentu (kondisional).
Tabel frekuensi parsial adalah tabel dua arah yang menunjukkan hubungan antara dua variabel (X dan Y) dengan mengendalikan (memfiksasi) nilai dari variabel ketiga (Z).
Tabel Frekuensi parsial untuk Z (Jenis Kelamin) Laki-Laki:* \[\begin{array}{|c|c|c|c|} \hline & \textbf{Seks Bebas (Ya)} & \textbf{Seks Bebas (Tidak)} & \textbf{Total} \\ \hline \textbf{Alkohol} & 130 & 100 & 230 \\ \textbf{Tidak Alkohol} & 30 & 90 & 120 \\ \hline \textbf{Total} & 160 & 190 & 350 \\ \hline \end{array}\]
Tabel Frekuensi parsial untuk Z (Jenis Kelamin) Perempuan:
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Seks Bebas (Ya)} & \textbf{Seks Bebas (Tidak)} & \textbf{Total} \\ \hline \textbf{Alkohol} & 70 & 90 & 160 \\ \textbf{Tidak Alkohol} & 20 & 90 & 110 \\ \hline \textbf{Total} & 90 & 180 & 270 \\ \hline \end{array}\]
Dengan menggunakan syntax R:
data3 <- array(c(130, 30, 100, 90, 70, 20, 90, 90 ),
dim = c(2, 2, 2),
dimnames = list(
Alkohol = c("Ya", "Tidak"),
Seks_Bebas = c("Ya", "Tidak"),
Jenis_Kelamin = c("Laki-Laki", "Perempuan")
))
# Ekstrak tabel parsial berdasarkan jenis kelamin
freq_parsial_laki <- data3[, , "Laki-Laki"]
freq_parsial_perempuan <- data3[, , "Perempuan"]
# Tampilkan hasil
freq_parsial_laki
## Seks_Bebas
## Alkohol Ya Tidak
## Ya 130 100
## Tidak 30 90
freq_parsial_perempuan
## Seks_Bebas
## Alkohol Ya Tidak
## Ya 70 90
## Tidak 20 90
Tabel frekuensi parsial adalah alat penting dalam analisis data kategorikal 3 arah, karena memungkinkan kita memahami dinamika hubungan antar variabel secara lebih terperinci dan kontekstual. Ini membantu menghindari generalisasi yang menyesatkan dari tabel total (agregat) dan mengarahkan pada interpretasi yang lebih akurat.
Tabel frekuensi marginal adalah tabel yang menunjukkan jumlah total dari suatu variabel dengan mengabaikan (mengakumulasi) variabel lain dalam tabel kontingensi. Pada tabel 3 arah yang melibatkan tiga variabel kategorikal (misalnya X, Y, dan Z), frekuensi marginal diperoleh dengan menjumlahkan frekuensi atas satu atau dua variabel.
Tabel Frekuensi Marginal untuk X (Alkohol) dan Y (Seks Bebas): \[\begin{array}{|c|c|c|c|} \hline & \textbf{Seks Bebas} & \textbf{Tidak Seks Bebas} & \textbf{Total} \\ \hline \textbf{Alkohol} & 200 & 190 & 390 \\ \textbf{Tidak Alkohol} & 50 & 180 & 230 \\ \hline \textbf{Total} & 250 & 370 & 620 \\ \hline \end{array}\]
Tabel Frekuensi Marginal untuk Z (Jenis Kelamin) dan Y (Seks Bebas):
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Seks Bebas} & \textbf{Tidak Seks Bebas} & \textbf{Total} \\ \hline \textbf{Laki-Laki} & 160 & 190 & 350 \\ \textbf{Perempuan} & 90 & 180 & 270 \\ \hline \textbf{Total} & 350 & 270 & 620 \\ \hline \end{array}\]
Dengan menggunakan syntax R:
data4 <- array(c(130, 30, 100, 90, 70, 20, 90, 90 ),
dim = c(2, 2, 2),
dimnames = list(
Alkohol = c("Ya", "Tidak"),
Seks_Bebas = c("Ya", "Tidak"),
Jenis_Kelamin = c("Laki-Laki", "Perempuan")
))
# Ekstrak tabel marginal
freq_marginal_X <- apply(data4,1,sum)
freq_marginal_Z <- apply(data4,3,sum)
freq_marginal_X
## Ya Tidak
## 390 230
freq_marginal_Z
## Laki-Laki Perempuan
## 350 270
Kesimpulan:
Tabel frekuensi marginal sangat penting dalam analisis data kategorikal karena memberikan informasi global tentang distribusi dua variabel. Memberikan gambaran umum distribusi dua variabel tanpa mempertimbangkan variabel ketiga. Membantu mengevaluasi pengaruh keseluruhan dari satu variabel terhadap yang lain.
Peluang bersama menunjukkan kemungkinan terjadinya ketiga peristiwa secara bersamaan:
\[ P(Z,X,Y) = \frac{{f(Z,X,Y)}}{N} \]
Sebagai contoh digunakan data pada bagian sebelumnya mengenai hubungan jenis kelamin, konsumsi alkohol, dan seks bebas. Peluang laki-laki, konumsi alkohol, dan perilaku seks bebas:
\[ P(Laki-Laki, Alkohol,Seks Bebas) = \frac{130}{620} = 0.20967 \]
Dengan menggunakan syntax R:
data <- data.frame(
Jenis_Kelamin = c(rep("Laki-laki", 2), rep("Perempuan", 2)),
Konsumsi_Alkohol = rep(c("Ya", "Tidak"), 2),
Seks_Bebas_Ya = c(130, 30, 70, 20),
Seks_Bebas_Tidak = c(100, 90, 90, 90),
Total = c(230, 120, 160, 110)
)
data$P_ZXY <- data$Seks_Bebas_Ya/sum(data$Total)
data
Tabel Probabilitas Bersama diperoleh dengan membagi semua nilai di dalam sel dengan jumlah keseluruhan sampel.
Tabel Probabilitas Bersama untuk Y (Seks Bebas) dan Tidak Seks Bebas:
\[\begin{array}{|c|c|cc|} \hline \textbf{Z (Jenis Kelamin)}& \textbf{X (Alkohol)} & \textbf{P(Y} & \textbf{X,Z)}\\ \hline \textbf{Laki-Laki} & Ya & 0.20967 & 0.16129 \\ \textbf{Laki-Laki} & Tidak & 0.04839 & 0.14516 \\ \textbf{Perempuan} & Ya & 0.11290 & 0.03226 \\ \textbf{Perempuan} & Tidak & 0.14516 & 0.14516 \\ \hline \end{array}\]
Dengan menggunakan syntax R:
data5 <- array(c(130, 30, 100, 90, 70, 20, 90, 90 ),
dim = c(2, 2, 2),
dimnames = list(
Alkohol = c("Ya", "Tidak"),
Seks_Bebas = c("Ya", "Tidak"),
Jenis_Kelamin = c("Laki-Laki", "Perempuan")
))
##Hitung probabilitas bersama
total <- sum(data5)
p_bersama <- data5/total
ftable(p_bersama)
## Jenis_Kelamin Laki-Laki Perempuan
## Alkohol Seks_Bebas
## Ya Ya 0.20967742 0.11290323
## Tidak 0.16129032 0.14516129
## Tidak Ya 0.04838710 0.03225806
## Tidak 0.14516129 0.14516129
Peluang marginal menunjukkan peluang terjadinya satu atau dua peristiwa, tanpa mempertimbangkan variabel lainnya. Untuk menghitung peluang marginal dapat dengan menjumlahkan probabilitas bersama.
\[ P(\text{Seks Bebas} = \text{Ya}) = \frac{130+30+70+20}{620} = \frac{250}{620} \approx 0.403225 \]
\[ P(\text{Seks Bebas} = \text{Tidak}) = \frac{100+90+90+90}{620} = \frac{370}{620} \approx 0.59677 \]
Dengan menggunakan syntax R:
P_Y1 <- sum(data$Seks_Bebas_Ya)/sum(data$Total)
P_Y1
## [1] 0.4032258
P_Y2 <- sum(data$Seks_Bebas_Tidak)/sum(data$Total)
P_Y2
## [1] 0.5967742
Tabel Peluang Marginal dihitung dengan menjumlahkan probabilitas bersama untuk setiap variabel yang tersisa setelah mengabaikan variabel lain.
Tabel marginal X (Alkohol) dan Y (Seks Bebas) \[\begin{array}{|c|c|cc|} \hline \textbf{X (Alkohol)} & \textbf{P(Y)} & \textbf{P(Tidak Y)}\\ \hline \textbf{Ya} & 0.30645 & 0.32258 \\ \textbf{Tidak} & 0.08065 & 0.29032\\ \hline \end{array}\]
\[\begin{array}{|c|c|cc|} \hline \textbf{Z (Jenis Kelamin)}& \textbf{P(Y)} & \textbf{P(Tidak Y)}\\ \hline \textbf{Laki-Laki} & 0.25806 & 0.30645 \\ \textbf{Perempuan} & 0.14516 & 0.29032 \\ \hline \end{array}\]
Dengan menggunakan syntax R:
data6 <- array(c(130, 30, 100, 90, 70, 20, 90, 90 ),
dim = c(2, 2, 2),
dimnames = list(
Alkohol = c("Ya", "Tidak"),
Seks_Bebas = c("Ya", "Tidak"),
Jenis_Kelamin = c("Laki-Laki", "Perempuan")
))
# Hitung probabilitas Bersama
total <- sum(data6)
p_bersama <- data6/total
# Hitung probabilitas marginal
marginal_X <- apply(p_bersama,1,sum)
marginal_Z <- apply(p_bersama,3,sum)
marginal_X
## Ya Tidak
## 0.6290323 0.3709677
marginal_Z
## Laki-Laki Perempuan
## 0.5645161 0.4354839
Kesimpulan
Tabel peluang marginal membantu memahami distribusi probabilitas dari masing-masing variabel tanpa mempertimbangkan variabel lainnya. Ini memberikan wawasan penting dalam analisis data kategori dengan tabel kontingensi tiga arah.
Peluang bersyarat didefinisikan sebagai: \[ P(Z|X,Y) = \frac{P(Z,X,Y)}{P(X,Y)} \]
\[\begin{align*} &\text{Total populasi: } N = 620 \\ \\ &\text{Probabilitas bersama: } \\ &P(Z = Laki-Laki, X = \text{Ya}, Y = \text{Ya}) = \frac{130}{620} \\ &P(X = \text{Ya}, Y = \text{Ya}) = \frac{130+70}{620} = \frac{200}{620} \\ \\ &\text{Peluang bersyarat:} \\ &P(Z = Laki-Laki \mid X = \text{Ya}, Y = \text{Ya}) = \frac{P(Z = Laki-Laki, X = \text{Ya}, Y = \text{Ya})}{P(X = \text{Ya}, Y = \text{Ya})} \\ &= \frac{130 / 620}{200 / 620} = \frac{130}{200} \approx 0.65 \\ \end{align*}\]
Dengan menggunakan syntax R:
P_X_Y1 <- sum(data$Seks_Bebas_Ya [data$Konsumsi_Alkohol == "Ya"]) / sum(data$Total)
P_Laki_given_X_Y1 <- (130 / 620) / P_X_Y1
P_Laki_given_X_Y1
## [1] 0.65
Kesimpulan:
Peluang bersama dihitung dengan membagi frekuensi dengan total populasi data
Peluang marginal dihitung dengan menjumlahkan probabilitas bersama
Peluang bersyarat dihitung dengan menggunakan rumus Bayes
Ukuran Asosiasi digunakan untuk mengukur kekuatan hubungan antara dua variabel kategori. Terdapat 3 ukuran asosiasi yang umum digunakan, yaitu:
Risk Difference
\[ RD = P(Y \mid X_{1},Z) - P(Y \mid X_{2},Z) \]
Relative Risk
\[ RR = \frac{P(Y \mid X_{1},Z)}{P(Y \mid X_{2},Z)} \]
\[ OR = \frac{P(Y \mid X_{1},Z)/(1-P(Y \mid X_{1},Z))}{P(Y \mid X_{2},Z)-(1-P(Y \mid X_{2},Z))} \]
Contoh Kasus:
Data menampilkan tabel kontingensi 3 arah, yang memperlihatkan hubungan pengaruh jenis kelamin, konsumsi minuman beralkohol terhadap perilaku seks bebas pada remaja. Data berjumlah 620 dengan proporsi 350 laki-laki, dan 270 perempuan.
\[\begin{array}{|c|c|c|c|} \hline \textbf{Jenis Kelamin} & \textbf{Alkohol} & \textbf{Seks Bebas: Ya} & \textbf{Seks Bebas: Tidak} \\ \hline Laki-Laki & Ya & 130 & 100 \\ Laki-Laki & Tidak & 30 & 90 \\ Perempuan & Ya & 70 & 90 \\ Perempuan & Tidak & 20 & 90 \\ \hline \end{array}\]
Tanpa memperthatikan variabel Z diperoleh tabel marginal sebagai berikut:
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Seks Bebas} & \textbf{Tidak Seks Bebas} & \textbf{Total} \\ \hline \textbf{Alkohol} & 200 & 190 & 390 \\ \textbf{Tidak Alkohol} & 50 & 180 & 230 \\ \hline \textbf{Total} & 250 & 370 & 620 \\ \hline \end{array}\]
Lakukan perhitungan manual:
Risk Difference:
\[ RD = 0.51282 - 0.21739 = 0.29174 \]
Relative Risk:
\[ RR = \frac{0.51282}{0.21739} = 2.35897 \]
Odds Ratio:
\[ OR = \frac{200 \times 180}{190 \times 50} = 3.78947 \]
Dengan menggunakan syntax R:
# Data untuk laki-laki
a_laki <- 130 # Alkohol, Seks Bebas
b_laki <- 100 # Alkohol, Tidak Seks Bebas
c_laki <- 30 # Tidak Alkohol, Seks Bebas
d_laki <- 90 # Tidak Alkohol, Tidak Seks Bebas
# Menghitung proporsi risiko
P_ekspos_laki <- a_laki / (a_laki + b_laki)
P_tidak_ekspos_laki <- c_laki / (c_laki + d_laki)
# Menghitung ukuran asosiasi
RD_laki <- P_ekspos_laki - P_tidak_ekspos_laki
RR_laki <- P_ekspos_laki / P_tidak_ekspos_laki
OR_laki <- (a_laki * d_laki) / (b_laki * c_laki)
# Menampilkan hasil
cat("Beda Risiko (RD) Laki-laki: ", RD_laki, "\n")
## Beda Risiko (RD) Laki-laki: 0.3152174
cat("Risiko Relatif (RR) Laki-laki: ", RR_laki, "\n")
## Risiko Relatif (RR) Laki-laki: 2.26087
cat("Odds Ratio (OR) Laki-laki: ", OR_laki, "\n")
## Odds Ratio (OR) Laki-laki: 3.9
Ketiga ukuran asosiasi ini menggambarkan bahwa terdapat hubungan positif antara konsumsi alkohol dan perilaku seks bebas. Namun, perlu diperhatikan lagi bahwa hubungan ini tidak selalu bersifat kasual, dan faktor lain juga dapat berkontribusi terhadap perilaku seks bebas.
Jika, nilai Z diperhatikan, maka akan didapatkan ukuran asosiasi sebagai berikut:
Risk Difference:
\[ RD = \frac{a}{a+b} - \frac{c}{c+d} \]
Relative Risk
\[ RR = \frac{\frac{a}{a+b}}{\frac{c}{c+d}} \]
Odds Ratio
\[ OR = \frac{a \times d}{b \times c} \]
Misalkan kita akan menghitung ukuran asosiasi untuk orang mengonsumsi alkohol dibandingkan dengan yang tidak mengonsumsi alkohol terhadap perliaku seks bebas (untuk jenis kelamin laki-laki).
\[\begin{array}{|c|c|c|c|} \hline \textbf{Jenis Kelamin} & \textbf{Alkohol} & \textbf{Seks Bebas: Ya} & \textbf{Seks Bebas: Tidak} \\ \hline Laki-Laki & Ya & 130 & 100 \\ Laki-Laki & Tidak & 30 & 90 \\ Perempuan & Ya & 70 & 90 \\ Perempuan & Tidak & 20 & 90 \\ \hline \end{array}\]
Perhitungan manual:
Kelompok Laki-Laki
Risk Difference :\[ RD = \frac{130}{130+100} - \frac{30}{30+90} = \frac{130}{230} - \frac{30}{120} = 0.3152174 \]
Relative Risk :
\[ RR = \frac{\frac{130}{130+100}}{\frac{30}{30+90}} = 2.26087 \]
\[ OR = \frac{130 \times 90}{100 \times 30} = 3.9 \]
Dengan menggunakan syntax R:
# Data untuk laki-laki
a_laki <- 130 # Alkohol, Seks Bebas
b_laki <- 100 # Alkohol, Tidak Seks Bebas
c_laki <- 30 # Tidak Alkohol, Seks Bebas
d_laki <- 90 # Tidak Alkohol, Tidak Seks Bebas
# Menghitung proporsi risiko
P_ekspos_laki <- a_laki / (a_laki + b_laki)
P_tidak_ekspos_laki <- c_laki / (c_laki + d_laki)
# Menghitung ukuran asosiasi
RD_laki <- P_ekspos_laki - P_tidak_ekspos_laki
RR_laki <- P_ekspos_laki / P_tidak_ekspos_laki
OR_laki <- (a_laki * d_laki) / (b_laki * c_laki)
# Menampilkan hasil
cat("Beda Risiko (RD) Laki-laki: ", RD_laki, "\n")
## Beda Risiko (RD) Laki-laki: 0.3152174
cat("Risiko Relatif (RR) Laki-laki: ", RR_laki, "\n")
## Risiko Relatif (RR) Laki-laki: 2.26087
cat("Odds Ratio (OR) Laki-laki: ", OR_laki, "\n")
## Odds Ratio (OR) Laki-laki: 3.9
Interpretasi:
Laki-laki yang mengonsumsi alkohol memiliki peluang 0,3152174 lebih tinggi untuk melakukan seks bebas dibandingkan laki-laki yang tidak mengonsumsi alkohol.
Laki-laki yang mengonsumsi alkohol 2,26087 kali lebih mungkin melakukan seks bebas dibandingkan laki-laki yang tidak mengonsumsi alkohol.
Laki-laki yang mengonsumsi alkohol 3,9 kali lebih mungkin melakukan seks bebas dibandingkan laki-laki yang tidak mengonsumsi alkohol.
Kelompok Perempuan
Risk Difference :\[ RD = \frac{70}{70+90} - \frac{20}{20+90} = \frac{70}{160} - \frac{20}{110} = 0.2556818 \]
Relative Risk :
\[ RR = \frac{\frac{70}{70+90}}{\frac{20}{20+90}} = 2.40625 \]
Odds Ratio :
\[ OR = \frac{70 \times 90}{20 \times 90} = 3.5 \]
Dengan menggunakan syntax R:
# Data untuk perempuan
a_perempuan <- 70 # Alkohol, Seks Bebas
b_perempuan <- 90 # Alkohol, Tidak Seks Bebas
c_perempuan <- 20 # Tidak Alkohol, Seks Bebas
d_perempuan <- 90 # Tidak Alkohol, Tidak Seks Bebas
# Menghitung proporsi risiko
P_ekspos_perempuan <- a_perempuan / (a_perempuan + b_perempuan)
P_tidak_ekspos_perempuan <- c_perempuan / (c_perempuan + d_perempuan)
# Menghitung ukuran asosiasi
RD_perempuan <- P_ekspos_perempuan - P_tidak_ekspos_perempuan
RR_perempuan <- P_ekspos_perempuan / P_tidak_ekspos_perempuan
OR_perempuan <- (a_perempuan * d_perempuan) / (b_perempuan * c_perempuan)
# Menampilkan hasil
cat("Beda Risiko (RD) Perempuan: ", RD_perempuan, "\n")
## Beda Risiko (RD) Perempuan: 0.2556818
cat("Risiko Relatif (RR) Perempuan: ", RR_perempuan, "\n")
## Risiko Relatif (RR) Perempuan: 2.40625
cat("Odds Ratio (OR) Perempuan: ", OR_perempuan, "\n")
## Odds Ratio (OR) Perempuan: 3.5
Interpretasi:
Perempuan yang mengonsumsi alkohol memiliki peluang 0,2556818 lebih tinggi untuk melakukan seks bebas dibandingkan perempuan yang tidak mengonsumsi alkohol.
Perempuan yang mengonsumsi alkohol 2,40625 kali lebih mungkin melakukan seks bebas dibandingkan perempuan yang tidak mengonsumsi alkohol.
Perempuan yang mengonsumsi alkohol 3,5 kali lebih mungkin melakukan seks bebas dibandingkan perempuan yang tidak mengonsumsi alkohol.
Terjadi ketika dua variabel menjadi independen setelah dikendalikan oleh variabel ketiga. Secara matematis, dua variabel X dan Y dikatakan independen secara kondisional terhadap variabel Z jika:
\[ P(X, Y \mid Z) = P(X \mid Z) P(Y \mid Z) \]
atau dalam bentuk frekuensi:
\[ \frac{n_{ijk}}{n_{k++}}=\frac{n_{i+k}}{n_{k++}} \times \frac{n_{+jk}}{n_{k++}} \]
Salah satu untuk melihat kondisional independen adalah uji chi-square dengan hipotesis:
Dengan statistik uji:
\[ \chi^2 = \sum \frac{(O_{ij} - E_{ij})^2}{E_{ij}} \]
di mana: - \(O_{ij}\) = Frekuensi yang diamati dalam kategori (observed frequency). - \(E_{ij}\) = Frekuensi yang diharapkan dalam kategori (expected frequency), dihitung sebagai:
\[ E_{ij} = \frac{(\text{Total baris}) \times (\text{Total kolom})}{\text{Total keseluruhan}} \]
Derajat kebebasan (df) dihitung dengan:
\[ df = (r-1)(c-1) \]
di mana:
\(r\) = jumlah baris (kategori konsumsi alkohol).
\(c\) = jumlah kolom (kategori seks bebas).
Melakukan Perhitungan:
Tabel kontingensi untuk laki-laki:
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Seks Bebas (Ya)} & \textbf{Seks Bebas (Tidak)} & \textbf{Total} \\ \hline \textbf{Alkohol} & 130 & 100 & 230 \\ \textbf{Tidak Alkohol} & 30 & 90 & 120 \\ \hline \textbf{Total} & 160 & 190 & 350 \\ \hline \end{array}\]
Hitung Expected Frequency (\(E\)):
\[ E_{11} = \frac{(230 \times 160)}{350} = 105.14 \]
\[ E_{12} = \frac{(230 \times 190)}{350} = 124.86 \]
\[ E_{21} = \frac{(120 \times 160)}{350} = 54.86 \]
\[ E_{22} = \frac{(120 \times 190)}{350} = 65.14 \]
Hitung \(\chi^2\):
\[ \chi^2 = \sum \frac{(O - E)^2}{E} \]
\[ \chi^2 = \frac{(130 - 105.14)^2}{105.14} + \frac{(100 - 124.86)^2}{124.86} + \frac{(30 - 54.86)^2}{54.86} + \frac{(90 - 65.14)^2}{65.14} \]
\[ \chi^2 = \frac{(24.86)^2}{105.14} + \frac{(-24.86)^2}{124.86} + \frac{(-24.86)^2}{54.86} + \frac{(24.86)^2}{65.14} \]
\[ \chi^2 = \frac{617.06}{105.14} + \frac{617.06}{124.86} + \frac{617.06}{54.86} + \frac{617.06}{65.14} \]
\[ \chi^2 = 5.87 + 4.94 + 11.25 + 9.48 \]
\[ \chi^2 = 31.54 \]
Tabel kontingensi untuk perempuan:
\[\begin{array}{|c|c|c|c|} \hline & \textbf{Seks Bebas (Ya)} & \textbf{Seks Bebas (Tidak)} & \textbf{Total} \\ \hline \textbf{Alkohol} & 70 & 90 & 160 \\ \textbf{Tidak Alkohol} & 20 & 90 & 110 \\ \hline \textbf{Total} & 90 & 180 & 270 \\ \hline \end{array}\]
Hitung Expected Frequency (\(E\)):
\[ E_{11} = \frac{(160 \times 90)}{270} = 53.33 \]
\[ E_{12} = \frac{(160 \times 180)}{270} = 106.67 \]
\[ E_{21} = \frac{(110 \times 90)}{270} = 36.67 \]
\[ E_{22} = \frac{(110 \times 180)}{270} = 73.33 \]
Hitung \(\chi^2\):
\[ \chi^2 = \sum \frac{(O - E)^2}{E} \]
\[ \chi^2 = \frac{(70 - 53.33)^2}{53.33} + \frac{(90 - 106.67)^2}{106.67} + \frac{(20 - 36.67)^2}{36.67} + \frac{(90 - 73.33)^2}{73.33} \]
\[ \chi^2 = \frac{(16.67)^2}{53.33} + \frac{(-16.67)^2}{106.67} + \frac{(-16.67)^2}{36.67} + \frac{(16.67)^2}{73.33} \]
\[ \chi^2 = \frac{278.09}{53.33} + \frac{278.09}{106.67} + \frac{278.09}{36.67} + \frac{278.09}{73.33} \]
\[ \chi^2 = 5.21 + 2.61 + 7.58 + 3.79 \]
\[ \chi^2 = 19.19 \]
Dengan menggunakan syntax R:
data_laki <- matrix(c(130, 100, 30, 90), nrow = 2, byrow = TRUE,
dimnames = list(c("Alkohol_Yes", "Alkohol_No"), c("Seks_Bebas_Yes", "Seks_Bebas_No")))
data_perempuan <- matrix(c(70, 90, 20, 90), nrow = 2, byrow = TRUE,
dimnames = list(c("Alkohol_Yes", "Alkohol_No"), c("Seks_Bebas_Yes", "Seks_Bebas_No")))
# Uji Chi-Square untuk kelompok laki-laki
chi_laki <- chisq.test(data_laki, correct=FALSE)
# Menampilkan hasil uji Chi-Square
print("Hasil Uji Chi-Square untuk Laki-laki:")
## [1] "Hasil Uji Chi-Square untuk Laki-laki:"
print(chi_laki)
##
## Pearson's Chi-squared test
##
## data: data_laki
## X-squared = 31.574, df = 1, p-value = 1.92e-08
# Uji Chi-Square untuk kelompok perempuan
chi_perempuan <- chisq.test(data_perempuan, correct=FALSE)
# Menampilkan hasil uji Chi-Square
print("Hasil Uji Chi-Square untuk Perempuan:")
## [1] "Hasil Uji Chi-Square untuk Perempuan:"
print(chi_perempuan)
##
## Pearson's Chi-squared test
##
## data: data_perempuan
## X-squared = 19.176, df = 1, p-value = 1.192e-05
Interpretasi:
Dalam analisis data kategorik, sering kali kita dihadapkan pada situasi di mana terdapat dua variabel utama yang ingin dianalisis hubungannya (misalnya variabel X dan Y), tetapi hubungan tersebut dapat dipengaruhi oleh variabel ketiga yang bersifat sebagai confounding atau variabel kontrol (Z). Untuk menangani situasi seperti ini, digunakanlah tabel kontingensi tiga arah, yang memungkinkan kita mengevaluasi asosiasi antara dua variabel sembari mengontrol efek dari variabel ketiga.
Sebagai contoh, kita ingin mengetahui apakah terdapat hubungan antara konsumsi alkohol (X) dengan kejadian perilau seks bebas (Y), dengan mempertimbangkan jenis kelamin sebagai variabel pengganggu (Z). Pendekatan yang digunakan adalah dengan membagi data ke dalam tabel parsial (yaitu tabel 2 × 2) untuk setiap strata dari variabel Z (misalnya, kelompok usia tertentu), serta membuat sebuah tabel marginal yang mengabaikan variabel Z.
Untuk menghitung kekuatan di tiap tabel parsial akan digunakan ukuran asosiasi yang umum digunakan, yakni Odds Ratio: \[ OR = \frac{a \cdot d}{b \cdot c} \] dengan:
Jika odds ratio parsial relatif konstan, kita dapat menghitung odds ratio bersama menggunakan estimasi Mantel-Haenszel.
\[ \text{OR}_{MH} = \frac{ \sum_{i=1}^{K} \frac{a_i d_i}{n_i} }{ \sum_{i=1}^{K} \frac{b_i c_i}{n_i} } \] Metode Mantel-Haenszel juga dapat digunakan jika ingin menguji apakah terdapat asosiasi antara variabel X dan Y dalam semua strata variabel Z.
Dalam analisis data kategorik, Independensi Bersyarat merupakan bagian penting yang memungkinkan kita untuk memahami hubungan antara dua variabel utama dengan mempertimbangkan adanya variabel ketiga yang berpotensi mempengaruhi hubungan tersebut. Secara konseptual, independensi bersyarat terjadi ketika dua variabel X dan Y tidak saling bergantung satu sama lain pada setiap tingkat strata variabel ketiga (Z). Secara matematis dapat dituliskan sebagai berikut:
\[ OR(X,Y \mid Z) = 1 \]
Artinya Odds Ratio antara X dan Y adalah 1 pada setiap strata Z, yang berarti tidak adanya hubungan atau ketergantungan antara keduanya dalam konteks strata tersebut.
Sebagai contoh misalkan seorang peneliti ingin mengetahui apakah hubungan variabel X (konsumsi alkohol) dan Y (seks bebas) tetap signifikan dengan mempertimbangkan variabel Z (jenis kelamin) sebagai variabel kontrol. Jika OR =1 atau mendekati 1 pada setiap strata Z maka independensi bersyarat berlaku dan variabel X dan Y independensi bersyarat terhadap variabel Z (jenis kelamin), dan jika nilai OR bervariasi secara nyata pada setiap strata Z maka independensi bersyarat tidak berlaku.
Pengujian Statistik Independensi Besyarat dilakukan dengan uji Cochran-Mantel-Haenszel dengan hipotesis sebagai berikut:
Hipotesis:
Statistik Uji CMH:
\[ CMH = \frac{\left[ \sum_k (n_{11k} - \mu_{11k}) \right]^2}{\sum_k \text{var}(n_{11k})} \]
Keterangan:
\[ \mu_{11k} = E(n_{11k}) = \frac{n_{1\cdot k} \cdot n_{\cdot 1k}}{n_{\cdot \cdot k}} \]
\[ \text{var}(n_{11k}) = \frac{n_{1\cdot k} \cdot n_{2\cdot k} \cdot n_{\cdot 1k} \cdot n_{\cdot 2k}}{n_{\cdot\cdot k}^2(n_{\cdot\cdot k}-1)} \]
Keputusan Uji:
Statistik uji CMH mengikuti distribusi chi-square dengan derajat kebebasan 1.
Tolak \(H_{0}\) jika,
\[ CMH > \chi^2_{(1)} \quad \text{atau p-value} < \alpha \]
Contoh Kasus:
\[\begin{array}{|c|c|c|c|} \hline \textbf{Jenis Kelamin} & \textbf{Alkohol} & \textbf{Seks Bebas: Ya} & \textbf{Seks Bebas: Tidak} \\ \hline Laki-Laki & Ya & 130 & 100 \\ Laki-Laki & Tidak & 30 & 90 \\ Perempuan & Ya & 70 & 90 \\ Perempuan & Tidak & 20 & 90 \\ \hline \end{array}\]
Hitung nilai ekspektasi:
\[ \mu_{11k} = \frac{n_{1\cdot k} \cdot n_{\cdot 1k}}{n_{\cdot\cdot k}} \]
\[\begin{align*} \mu_{11(1)} &= \frac{230 \times 160}{350} = 105.143 \\ \mu_{11(2)} &= \frac{160 \times 90}{270} = 53.33 \\ \end{align*}\]
Hitung nilai varians:
\[ \text{Var}(n_{11k}) = \frac{n_{1\cdot k} \cdot n_{2\cdot k} \cdot n_{\cdot 1k} \cdot n_{\cdot 2k}}{n_{\cdot\cdot k}^2 (n_{\cdot\cdot k} - 1)} \]
\[\begin{align*} \text{Var}(n_{11(1)}) &= \frac{230 \times 120 \times 160 \times 190}{350^2 \times 349} = 19.626 \\ \text{Var}(n_{11(2)}) &= \frac{160 \times 110 \times 90 \times 180}{270^2 \times 269} = 14.539 \\ \end{align*}\]
Hitung statistik uji CMH:
\[ X^2_{CMH} = \frac{\left( \sum_k (n_{11k} - \mu_{11k}) \right)^2}{\sum_k \text{Var}(n_{11k})} \]
\[\begin{align*} X^2_{CMH} &= \frac{1724.227}{19.626 + 14.539} \\ &= \frac{1724.227}{34.165} \\ &= 50.468 \end{align*}\]
Nilai tersebut dibandingkan dengan distribusi \(\chi^2\) dengan derajat kebebasan 1.
Dengan \(\alpha = 0.05\), nilai kritis \(\chi^2_{1,0.05} = 3.841\). Karena \(174.087 > 3.841\), maka hipotesis nol ditolak, yang berarti ada hubungan signifikan antara konsumsi alkohol dan perilaku seks bebas setelah mengendalikan variabel jenis kelamin. Ini menunjukkan bahwa konsumsi alkohol dan perilaku seks bebas tidak bersifat conditional independence setelah jenis kelamin dikendalikan.
Dengan menggunakan syntax R:
# Membuat array 2x2x2 (2 lapisan untuk jenis kelamin)
data_cmh <- array(c(130, 30, 100, 90, # Laki-laki
70, 20, 90, 90), # Perempuan
dim = c(2, 2, 2),
dimnames = list(
"Alkohol" = c("Ya", "Tidak"),
"Seks Bebas" = c("Ya", "Tidak"),
"Jenis Kelamin" = c("Laki-laki", "Perempuan")
))
# Menampilkan tabel 2x2x2 untuk CMH
print("Tabel 2x2x2 untuk Uji CMH:")
## [1] "Tabel 2x2x2 untuk Uji CMH:"
print(data_cmh)
## , , Jenis Kelamin = Laki-laki
##
## Seks Bebas
## Alkohol Ya Tidak
## Ya 130 100
## Tidak 30 90
##
## , , Jenis Kelamin = Perempuan
##
## Seks Bebas
## Alkohol Ya Tidak
## Ya 70 90
## Tidak 20 90
# Melakukan Uji Cochran-Mantel-Haenszel (CMH)
library(vcd)
## Warning: package 'vcd' was built under R version 4.3.3
## Loading required package: grid
cmh_test <- mantelhaen.test(data_cmh, correct=FALSE)
# Menampilkan hasil uji CMH
print("Hasil Uji Cochran-Mantel-Haenszel:")
## [1] "Hasil Uji Cochran-Mantel-Haenszel:"
print(cmh_test)
##
## Mantel-Haenszel chi-squared test without continuity correction
##
## data: data_cmh
## Mantel-Haenszel X-squared = 50.468, df = 1, p-value = 1.211e-12
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 2.565755 5.408009
## sample estimates:
## common odds ratio
## 3.725
Didapatkan nilai p-value (1.211e-12) < 0.05 maka Tolak H0. Jadi, dapat disimpulkan bahwa terdapat hubungan signifikan antara konsumsi alkohol (X) dan seks bebas (Y) setelah mengontrol faktor jenis kelamin (Z).
Odds Ratio (OR) adalah ukuran asosiasi dalam tabel kontingensi 2x2 yang menunjukkan seberapa besar peluang suatu kejadian terjadi pada satu kelompok dibandingkan dengan kelompok lain. Namun, ketika data dikategorikan menjadi beberapa strata berdasarkan variabel perancu (seperti usia), kita tidak cukup hanya menghitung odds ratio dari satu tabel. Untuk itu, kita gunakan Odds Ratio Bersama atau Common Odds Ratio menggunakan pendekatan Mantel-Haenszel (MH).
\[ \hat{\theta}{MH} = \frac{\sum_{k=1}^{K} \left( \frac{n_{11k} n_{22k}}{n_{..k}} \right)}{\sum_{k=1}^{K} \left( \frac{n_{12k} n_{21k}}{n_{..k}} \right)} \]
Standar error log Odds Ratio bersama:
\[ \hat{\sigma}^{2}[\log(\hat{\theta}{MH})] = \frac{\sum (n{11k} + n_{12k})(n_{11k}n_{22k})/n_{..k}^{2}}{2(\sum n_{11k}n_{12k}/n_{..k})^{2}} + \frac{\sum(n_{11k} + n_{22k})(n_{11k} + n_{12k}) + (n_{12k} + n_{21k})(n_{11k} + n_{22k})]/n_{..k}^{2}}{2(\sum n_{11k}n_{12k}/n_{..k})^{2}} + \frac{\sum(n_{12k} + n_{21k})(n_{12k}n_{21k})/n_{..k}^{2}}{2(\sum n_{12k}n_{21k}/n_{..k})^{2}} \]
Interval Kepercayaan Odds Ratio Bersama:
\[ \log(\hat{\theta}{MH}) \pm Z{\alpha} \hat{\sigma}[\log(\hat{\theta}_{MH})] \]
Implementasi pada data konusmsi alkohol dan perilaku seks bebas:
\[ \begin{array}{|c|c|c|c|} \hline \textbf{Jenis Kelamin} & \textbf{Alkohol} & \textbf{Seks Bebas: Ya} & \textbf{Seks Bebas: Tidak} \\ \hline Laki-Laki & Ya & 130 & 100 \\ Laki-Laki & Tidak & 30 & 90 \\ Perempuan & Ya & 70 & 90 \\ Perempuan & Tidak & 20 & 90 \\ \hline \end{array} \]
Taksiran Odds Ratio Bersama:
\[\begin{align*} \text{Pembilang} &= \frac{(130)(90)}{350} + \frac{(70)(90)}{270} \\ &= 33.429 + 23.333 = 56.762 \\ \\ \text{Penyebut} &= \frac{(100)(30)}{350} + \frac{(90)(20)}{270} \\ &= 8.571 + 6.667 = 15.238 \\ \\ \hat{\theta}_{MH} &= \frac{56.762}{15.238} = 3.725 \end{align*}\]
Standar error log Odds Ratio:
Misalkan nilai standar error log Odds Ratio Bersama = 0.05.
Interval Kepercayaan Odds Ratio Bersama:
\[\begin{align*} \log(3.725) &\pm 1.96 \times 0.06 \\ 1.315 &\pm 0.098 \\ \Rightarrow \text{CI log} &= [1.217,\ 1.413] \\ \Rightarrow \text{CI OR} &= [e^{1.217},\ e^{1.413}] = [3.378,\ 4.108] \end{align*}\]
Interpretasi Hasil:
Odds ratio bersama menunjukkan bahwa konsumen alkohol 3,75 kali lebih mungkin melakukan seks bebas dibandingkan yang tidak mengonsumsi alkohol.
Standar error = 0.05 artinya variabilitas kecil pada estimasi cukup ratio
Interval Kepercayaan (3.378 hingga 4.108) menunjukkan bahwa kita cukup yakin bahwa odds ratio yang sebenarnya berada dalam interval ini.
Interval Kepercayaan tidak mencakup 1, hasil ini signifikan secara statistik, sekaligus mendukung adanya hubungan antara konsumsi alkohol dan perilaku seks bebas.
Dengan menggunakan syntax R:
# Membuat array 2x2x2 (2 lapisan untuk jenis kelamin)
library(epitools)
##
## Attaching package: 'epitools'
## The following object is masked from 'package:vcd':
##
## oddsratio
data_cmh <- array(c(130, 30, 100, 90, # Laki-laki
70, 20, 90, 90), # Perempuan
dim = c(2, 2, 2),
dimnames = list(
"Alkohol" = c("Ya", "Tidak"),
"Seks Bebas" = c("Ya", "Tidak"),
"Jenis Kelamin" = c("Laki-laki", "Perempuan")
))
# Menampilkan tabel 2x2x2 untuk CMH
print("Tabel 2x2x2 untuk Uji CMH:")
## [1] "Tabel 2x2x2 untuk Uji CMH:"
print(data_cmh)
## , , Jenis Kelamin = Laki-laki
##
## Seks Bebas
## Alkohol Ya Tidak
## Ya 130 100
## Tidak 30 90
##
## , , Jenis Kelamin = Perempuan
##
## Seks Bebas
## Alkohol Ya Tidak
## Ya 70 90
## Tidak 20 90
# Melakukan Uji Cochran-Mantel-Haenszel (CMH)
cmh_test <- mantelhaen.test(data_cmh, correct=FALSE)
# Menampilkan hasil uji CMH
print("Hasil Uji Cochran-Mantel-Haenszel:")
## [1] "Hasil Uji Cochran-Mantel-Haenszel:"
print(cmh_test)
##
## Mantel-Haenszel chi-squared test without continuity correction
##
## data: data_cmh
## Mantel-Haenszel X-squared = 50.468, df = 1, p-value = 1.211e-12
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 2.565755 5.408009
## sample estimates:
## common odds ratio
## 3.725
Dalam analisis tabel kontingensi 3 arah, kita ingin mengetahui apakah hubungan antara variabel X dan Y konsisten di seluruh strata dari variabel ketiga Z. Untuk mengetahui hal tersebut kita dapat menggunakan Uji Homogenitas Odds Ratio yang salah satu metode nya menggunakan statistik uji Breslow-Day.
Asosiasi dikatakan homogen apabila nilai odds ratio dari setiap strata memiliki nilai yang sama. Secara matematis kondisi ini dituliskan sebagai:
\[ \theta_{xy(1)} = \theta_{xy(2)} = \cdots = \theta_{xy(K)} \]
Jika kondisi ini terpenuhi maka dapat disimpulkan bahwa tidak terdapat interaksi atau hubungan X dan Y tidak bergantung pada variabel Z. Namun, jika kondisi tersebut tidak terpenuhi maka terindikasi bahwa hubungan X dan Y dipengaruhi oleh Z.
Untuk mengetahui hal tersebut akan dilakukan pengujian statistik menggunakan statistik uji Breslow-Day dengan hipotesis sebagai:
Hipotesis Nol (\(H_0\)):
Tidak ada perbedaan odds ratio antar strata, yaitu:
\[ \theta_{xy(1)} = \theta_{xy(2)} = \cdots = \theta_{xy(K)} \]
Hipotesis Alternatif (\(H_1\)):
Setidaknya terdapat satu strata di mana nilai odds ratio berbeda dari yang lain.
Statistik Uji Breslow-Day:
\[ X^2_{HBD} = \sum_{j=1}^{K} \frac{(a_j - \tilde{a}j)^2}{\widehat{\text{Var}}(a_j | \hat{OR}{MH})} \]
Statistik uji Breslow-Day mengikuti distribusi chi square dengan derajat kebebasan = K-1
Implementasi pada contoh kasus hubungan konsumsi alkohol dan perilaku seks bebas:
\[ \begin{array}{|c|c|c|c|} \hline \textbf{Jenis Kelamin} & \textbf{Alkohol} & \textbf{Seks Bebas: Ya} & \textbf{Seks Bebas: Tidak} \\ \hline Laki-Laki & Ya & 130 & 100 \\ Laki-Laki & Tidak & 30 & 90 \\ Perempuan & Ya & 70 & 90 \\ Perempuan & Tidak & 20 & 90 \\ \hline \end{array} \]
Hitung Odds Ratio setiap strata:
\[ OR = \frac{a \times d}{b \times c} \]
Untuk setiap kelompok jenis kelamin:
\[\begin{align*} OR_1 &= \frac{(130 \times 90)}{(100 \times 30)} = \frac{11700}{3000} = 3.9 \\ OR_2 &= \frac{(70 \times 90)}{(90 \times 20)} = \frac{6300}{1800} = 3.5 \\ \end{align*}\]
Hitung Odds Ratio Bersama:
\[\begin{align*} \text{Pembilang} &= \frac{(130)(90)}{350} + \frac{(70)(90)}{270} \\ &= 33.429 + 23.333 = 56.762 \\ \\ \text{Penyebut} &= \frac{(100)(30)}{350} + \frac{(90)(20)}{270} \\ &= 8.571 + 6.667 = 15.238 \\ \\ \hat{\theta}_{MH} &= \frac{56.762}{15.238} = 3.725 \end{align*}\]
Hitung Ekspektasi \(\bar{a}_j\):
\[ -m_{j1} n_{j1} \hat{OR}{MH} + (n{j2} - m_{j1} + \hat{OR}{MH}(m{j1} + m_{j2})) \bar{a}j + (1 - \hat{OR}{MH}) \bar{a}_j^2 = 0 \]
\[\begin{align*} \bar{a}_1 &= \frac{m_{11} \cdot n_{11} \cdot \theta}{\theta(n_{11} + m_{21}) + m_{11} + n_{21}} \\ &= \frac{230 \cdot 160 \cdot 3.72}{3.72(160 + 120) + 230 + 190} \\ &= \frac{136896}{1461.6} \approx 93.7 \end{align*}\]
\[\begin{align*} \bar{a}_2 &= \frac{m_{1,2} \cdot n_{1,2} \cdot \theta}{\theta(n_{1,2} + m_{2,2}) + m_{1,2} + n_{2,2}} \\ &= \frac{160 \cdot 90 \cdot 3.72}{3.72(90 + 110) + 160 + 180} \\ &= \frac{53568}{744 + 340} = \frac{53568}{1084} \approx 49.42 \end{align*}\]
Hitung varians: \[ \text{Var}(a_j | \hat{OR}_{MH}) = \left( \frac{1}{a_j} + \frac{1}{b_j} + \frac{1}{c_j} + \frac{1}{d_j} \right)^{-1} \]
\[\begin{align*} \mathrm{Var}(a_1 \mid \hat{OR}_{MH}) &= \left( \frac{1}{a_1} + \frac{1}{b_1} + \frac{1}{c_1} + \frac{1}{d_1} \right)^{-1} \\ &= \left( \frac{1}{130} + \frac{1}{100} + \frac{1}{30} + \frac{1}{90} \right)^{-1} \\ &= \left( 0.00769 + 0.01 + 0.03333 + 0.01111 \right)^{-1} \\ &= \left( 0.06213 \right)^{-1} \approx 16.10 \end{align*}\]
\[\begin{align*} \mathrm{Var}(a_2 \mid \hat{OR}_{MH}) &= \left( \frac{1}{a_2} + \frac{1}{b_2} + \frac{1}{c_2} + \frac{1}{d_2} \right)^{-1} \\ &= \left( \frac{1}{70} + \frac{1}{90} + \frac{1}{20} + \frac{1}{90} \right)^{-1} \\ &= \left( 0.01429 + 0.01111 + 0.05 + 0.01111 \right)^{-1} \\ &= \left( 0.08651 \right)^{-1} \approx 11.56 \end{align*}\]
Hitung Statistik Uji Breslow-Day:
\[ \chi^2_{BD} = \sum_j \frac{(a_j - \bar{a}_j)^2}{\mathrm{Var}(a_j \mid \hat{OR}_{MH})} \]
\[\begin{align*} \chi^2_{BD} &= \frac{(130 - 93.7)^2}{16.10} + \frac{(70 - 56.3)^2}{11.56} \\&= \frac{(36.3)^2}{16.10} + \frac{(13.7)^2}{11.56} \\&= \frac{1317.69}{16.10} + \frac{187.69}{11.56} \\&\approx 81.86 + 16.23 = \boxed{98.09} \end{align*}\]
Terapkan Koreksi Tarone:
\[ X^2_{HBDT} = X^2_{HBD} - \frac{\sum_j (a_j - \bar{a}_j)^2}{\sum_j \text{Var}(a_j)} \]
\[ X^2_{HBDT} = 98.09 - \frac{1505.38}{27.66} = 98.09 \times 54.41 = 43.68 \]
Hitung P-Value:
\[ p = 1 - \chi^2(X^2_{HBDT}, df = 1) \]
\[ p = 1 - \chi^2(43.68, df = 1) \approx 3.87e-11 \]
Dengan menggunakan syntax R:
library(DescTools)
## Warning: package 'DescTools' was built under R version 4.3.3
data_bdt <- array(c(130, 30, 100, 90, 70, 20, 90, 90),
dim = c(2, 2, 2),
dimnames = list(
"Alkohol" = c("Ya", "Tidak"),
"Seks Bebas" = c("Ya", "Tidak"),
"Jenis Kelamin" = c("Laki-laki", "Perempuan")
))
print(data_bdt)
## , , Jenis Kelamin = Laki-laki
##
## Seks Bebas
## Alkohol Ya Tidak
## Ya 130 100
## Tidak 30 90
##
## , , Jenis Kelamin = Perempuan
##
## Seks Bebas
## Alkohol Ya Tidak
## Ya 70 90
## Tidak 20 90
breslow_test <- BreslowDayTest(data_bdt)
print(breslow_test)
##
## Breslow-Day test on Homogeneity of Odds Ratios
##
## data: data_bdt
## X-squared = 0.078803, df = 1, p-value = 0.7789
Interpretasi:
Jika p-value < alpha (0.05), maka Tolak H0.
Jika p-value > alpha (0.05), maka Terima H0 artinya hubungan variabel X dan Y tidak dipengaruhi oleh variabel ketiga Z
Generalized Linear Model (GLM) adalah sebuah kerangka statistik yang digunakan untuk memodelkan hubungan antara variabel dependen dengan satu atau lebih variabel independen. GLM adalah perluasan dari model linier klasik yang memungkinkan untuk menggunakan berbagai jenis distribusi untuk variabel dependen, bukan hanya distribusi normal.
GLM terdiri dari tiga komponen utama:
Distribusi probabilitas untuk variabel
dependen
GLM memungkinkan penggunaan distribusi selain distribusi normal, seperti
distribusi binomial, Poisson, dan lainnya. Distribusi ini tergantung
pada sifat data yang dianalisis.
Fungsi link
Fungsi link menghubungkan nilai ekspektasi dari variabel dependen (\(E[Y]\)) dengan kombinasi linier dari
variabel independen. Fungsi link ini mengatur bagaimana hubungan antara
prediktor dan outcome.
Kombinasi linier prediktor
Kombinasi linier prediktor adalah bentuk linear dari variabel independen
yang dihubungkan dengan koefisien model. Ini adalah bagian dari model
yang memberikan prediksi dasar.
Exponential family adalah sekumpulan distribusi probabilitas yang memiliki bentuk fungsi kepadatan atau fungsi massa probabilitas yang dapat dituliskan dalam bentuk eksponensial. Keluarga distribusi ini mencakup banyak distribusi penting dalam statistika, termasuk distribusi normal, binomial, Poisson, dan banyak lainnya.
Distribusi yang termasuk dalam exponential family memiliki bentuk umum sebagai berikut:
\[ f(y; \theta, \phi) = \exp\left( \frac{y\theta - b(\theta)}{a(\phi)} + c(y, \phi) \right) \]
di mana:
\(\theta\) adalah natural parameter
\(\phi\) adalah dispersion parameter
\(a(\phi)\), \(b(\theta)\), dan \(c(y, \phi)\) adalah fungsi-fungsi tertentu yang berbeda untuk setiap distribusi
Distribusi yang termasuk dalam keluarga ini mencakup banyak distribusi umum seperti:
Distribusi Normal
Distribusi Binomial
Distribusi Poisson
Distribusi Gamma
Untuk menunjukkan bahwa suatu distribusi termasuk ke dalam exponential family, kita perlu menuliskan bentuk fungsi kepadatannya (PDF/PMF) dan menyusunnya agar sesuai dengan bentuk umum di atas.
Langkah-langkah:
Contoh: Distribusi Poisson
Distribusi Poisson memiliki fungsi probabilitas:
\[ f(y; \lambda) = \frac{\lambda^y e^{-\lambda}}{y!} \]
Manipulasi ke Bentuk Exponential Family
\[ f(y; \lambda) = \exp\left( y \log \lambda - \lambda - \log y! \right) \]
Kita bisa tuliskan dalam bentuk:
\[ f(y; \theta) = \exp\left( y\theta - e^{\theta} - \log y! \right) \]
Dengan substitusi \(\theta = \log \lambda\), maka \(\lambda = e^\theta\)
Identifikasi:
Kesimpulan: Distribusi Poisson adalah anggota dari exponential family.
Persamaan regresi logistik merupakan regresi linear, di mana nilai input dikombinasikan secara linear dengan sebuah bobot, yaitu menggunakan fungsi linear. Namun, regresi logistik membatasi hasil prediksi menjadi nilai antara 0 dan 1, dengan menggunakan fungsi aktivasi sigmoid.
Output yang diperoleh akan direpresentasikan dalam rentang nilai 0 hingga 1 dan sangatlah bersifat kurva S (S-shaped).
Regresi logistik menggambarkan hubungan antara satu atau lebih variabel independen dan menghasilkan nilai dalam dua kelas kategorik (biner). Model ini sangat cocok digunakan pada data dengan dua kelas keluaran (misalnya: ya/tidak, positif/negatif, sukses/gagal, dll).
Dengan kata lain, apabila di depan variabel bebas digunakan angka 1 mewakili kelas positif. Regresi logistik akan menghasilkan probabilitas untuk prediksi output kelas 1, dan variabel bebas tersebut akan dituliskan sebagai \(x = (x_1, x_2, ..., x_n)\).
Model secara umum berupaya memodelkan probabilitas bahwa \(Y = 1\) sebagai fungsi dari prediktor.
Beberapa contoh penerapan klasifikasi biner ini dapat ditemukan dalam berbagai kehidupan sehari-hari, seperti:
Keunggulan Utama Regresi Logistik
Regresi logistik lebih mudah dibandingkan dengan pendekatan neural network karena modelnya lebih sederhana dan interpretasinya lebih jelas. Beberapa keuntungan lainnya adalah:
Fungsi Aktivasi Sigmoid
Untuk mengubah output linear menjadi probabilitas, digunakan fungsi sigmoid. Fungsi sigmoid adalah:
\[ f(x) = \frac{1}{1 + e^{-x}} \]
Dengan fungsi ini, hasil prediksi akan bernilai antara 0 hingga 1. Beberapa interpretasi:
Fungsi sigmoid inilah yang menjadikan regresi logistik berbeda dari model linear biasa.Plot kurva sigmoid sebagai berikut:
# Simulasi data untuk regresi logistik - Peluang Lulus Ujian
set.seed(123)
n <- 120
jam_belajar <- seq(0, 10, length.out = n)
# Ubah parameter regresi
intercept <- -3
slope <- 0.8
# Hitung log odds dan probabilitas
log_odds <- intercept + slope * jam_belajar
prob_lulus <- 1 / (1 + exp(-log_odds))
# Simulasi label (lulus atau tidak)
status_lulus <- rbinom(n, 1, prob_lulus)
# Buat data frame
data2 <- data.frame(
jam_belajar = jam_belajar,
status_lulus = status_lulus,
prob_lulus = prob_lulus
)
# Visualisasi
plot(jam_belajar, status_lulus, pch = 16, col = "darkorange",
xlab = "Jam Belajar", ylab = "Status Lulus (0 = Gagal, 1 = Lulus)",
main = "Simulasi Regresi Logistik: Peluang Lulus Ujian")
lines(jam_belajar, prob_lulus, col = "darkgreen", lwd = 2)
abline(h = 0.5, col = "red", lty = 2)
legend("bottomright",
legend = c("Data Biner", "Kurva Logistik", "Ambang 0.5"),
col = c("darkorange", "darkgreen", "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 antara variabel prediktor dan probabilitas output. Pendekatan ini efektif untuk klasifkasi biner seperti deteksi penyakit, email spam, dan prediksi ya/tidak
Fungsi Link Logit dalam Regresi Logistik
Fungsi link (logit)
Fungsi link logit digunakan dalam regresi logistik untuk
menghubungkan nilai rata-rata dari respon (\(\mu\)) dengan kombinasi linear dari
prediktor. Fungsi ini didefinisikan sebagai:
\[ g(\mu) = \log \left( \frac{\mu}{1 - \mu} \right) \]
Fungsi ini merepresentasikan log-odds dari probabilitas keberhasilan \(\mu\). Artinya, ia mengubah nilai probabilitas yang berada pada rentang \((0, 1)\) menjadi nilai yang tidak dibatasi dalam rentang \((-\infty, +\infty)\).
Model Regresi Logistik
Model regresi logistik secara umum dituliskan sebagai:
\[ \log \left( \frac{\mu}{1 - \mu} \right) = X\beta \]
di mana: - \(\mu\) adalah probabilitas keberhasilan (misalnya, peluang sukses, lulus, positif, dll), - \(X\) adalah matriks variabel independen (fitur/prediktor), - \(\beta\) adalah vektor koefisien regresi.
Fungsi inverse link
Untuk memperoleh nilai probabilitas \(\mu\) dari kombinasi linear \(X\beta\), digunakan fungsi inverse dari logit, yaitu:
\[ \mu = \frac{\exp(X\beta)}{1 + \exp(X\beta)} \]
Fungsi ini dikenal sebagai fungsi sigmoid, yang memetakan nilai dari seluruh garis real ke interval \((0, 1)\), cocok untuk merepresentasikan probabilitas dalam regresi logistik.
Estimasi Parameter Dalam regresi logistik, tujuan utama dari proses estimasi adalah untuk menentukan parameter \(\beta\) yang memaksimalkan kesesuaian model terhadap data yang diamati. Karena model regresi logistik bukan linear dalam artian biasa, maka metode kuadrat terkecil (least squares) tidak dapat digunakan secara langsung.
Sebagai gantinya, digunakan metode Maximum Likelihood Estimation (MLE) atau Estimasi Likelihood Maksimum.
Prinsip dasar dari MLE adalah mencari nilai parameter \(\beta\) yang memaksimalkan fungsi likelihood, yaitu peluang bahwa model menghasilkan data seperti yang diamati.
Misalkan kita punya data pengamatan \((x_i, y_i)\), dengan: - \(y_i \in \{0, 1\}\) - \(\mu_i = P(y_i = 1 | x_i) = \frac{\exp(x_i^T \beta)}{1 + \exp(x_i^T \beta)}\)
Fungsi likelihood untuk seluruh data (dengan asumsi independen) dapat dituliskan sebagai:
\[ L(\beta) = \prod_{i=1}^{n} \mu_i^{y_i} (1 - \mu_i)^{1 - y_i} \]
Agar lebih mudah dihitung, digunakan log-likelihood:
\[ \ell(\beta) = \sum_{i=1}^{n} \left[ y_i \log(\mu_i) + (1 - y_i) \log(1 - \mu_i) \right] \]
Substitusi \(\mu_i = \frac{\exp(x_i^T \beta)}{1 + \exp(x_i^T \beta)}\), maka:
\[ \ell(\beta) = \sum_{i=1}^{n} \left[ y_i (x_i^T \beta) - \log(1 + \exp(x_i^T \beta)) \right] \]
Proses estimasi dilakukan dengan mengoptimalkan fungsi log-likelihood ini, biasanya dengan metode numerik, seperti:
# 1. Simulasi data (ulang agar sesuai)
set.seed(123)
n <- 120
jam_belajar <- seq(0, 10, length.out = n)
# Parameter logit model
intercept <- -3
slope <- 0.8
# Hitung probabilitas lulus
log_odds <- intercept + slope * jam_belajar
prob_lulus <- 1 / (1 + exp(-log_odds))
# Simulasi label biner (0 = gagal, 1 = lulus)
status_lulus <- rbinom(n, 1, prob_lulus)
# Buat data frame
data_logistik <- data.frame(
jam_belajar = jam_belajar,
status_lulus = status_lulus
)
# 2. Estimasi model regresi logistik
model_logistik <- glm(status_lulus ~ jam_belajar, data = data_logistik, family = binomial)
# 3. Lihat hasil estimasi parameter
summary(model_logistik)
##
## Call:
## glm(formula = status_lulus ~ jam_belajar, family = binomial,
## data = data_logistik)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.5339 0.5569 -4.550 5.35e-06 ***
## jam_belajar 0.7458 0.1326 5.626 1.85e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 155.387 on 119 degrees of freedom
## Residual deviance: 92.544 on 118 degrees of freedom
## AIC: 96.544
##
## Number of Fisher Scoring iterations: 5
data_logistik$pred_prob <- predict(model_logistik, type = "response")
data_logistik$pred_logit <- predict(model_logistik, type = "link")
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.3.3
# --- 3. Plot Kurva Sigmoid (Probabilitas Lulus) ---
ggplot(data_logistik, aes(x = jam_belajar, y = status_lulus)) +
geom_point(color = "gray60", alpha = 0.6, size = 2) +
geom_line(aes(y = pred_prob), color = "blue", size = 1.2) +
geom_hline(yintercept = 0.5, linetype = "dashed", color = "red") +
labs(title = "Kurva Sigmoid Regresi Logistik",
x = "Jam Belajar",
y = "Probabilitas / Status Lulus") +
theme_minimal()
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
Evaluasi Model
# --- 1. Prediksi klasifikasi dengan threshold 0.5 ---
data_logistik$pred_class <- ifelse(data_logistik$pred_prob > 0.5, 1, 0)
# Confusion Matrix
table_pred <- table(Prediksi = data_logistik$pred_class, Aktual = data_logistik$status_lulus)
print(table_pred)
## Aktual
## Prediksi 0 1
## 0 31 10
## 1 11 68
# Hitung akurasi
akurasi <- sum(diag(table_pred)) / sum(table_pred)
cat("Akurasi model:", round(akurasi, 4), "\n")
## Akurasi model: 0.825
# --- 2. Nilai AIC dan Deviance ---
cat("AIC:", AIC(model_logistik), "\n")
## AIC: 96.5443
cat("Null Deviance:", model_logistik$null.deviance, "\n")
## Null Deviance: 155.3872
cat("Residual Deviance:", model_logistik$deviance, "\n")
## Residual Deviance: 92.5443
# --- 3. Pseudo R-squared (McFadden's R²) ---
pseudoR2 <- 1 - (model_logistik$deviance / model_logistik$null.deviance)
cat("Pseudo R² (McFadden):", round(pseudoR2, 4), "\n")
## Pseudo R² (McFadden): 0.4044
Regresi Poisson adalah salah satu model dalam keluarga Generalized Linear Model (GLM) yang digunakan untuk memodelkan data diskrit berupa hitungan (jumlah kejadian dalam suatu interval waktu, ruang, atau unit lain).
Contoh kasus:
Jumlah pelanggan yang datang ke toko dalam 1 hari
Jumlah kecelakaan di suatu jalan dalam sebulan
Jumlah email spam yang diterima per minggu
Fungsi link: \[ \log(\lambda) = \beta_0 + \beta_1x \]
Simulasi Data dan Evaluasi Model:
set.seed(1)
n <- 100
x <- rnorm(n, mean = 0, sd = 1)
lambda <- exp(0.5 + 0.4 * x)
y <- rpois(n, lambda)
data_poisson <- data.frame(x = x, y = y)
model_poisson <- glm(y ~ x, data = data_poisson, family = poisson)
summary(model_poisson)
##
## Call:
## glm(formula = y ~ x, family = poisson, data = data_poisson)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.30986 0.09134 3.392 0.000693 ***
## x 0.50125 0.09009 5.564 2.64e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 137.41 on 99 degrees of freedom
## Residual deviance: 106.15 on 98 degrees of freedom
## AIC: 303.81
##
## Number of Fisher Scoring iterations: 5
exp(coef(model_poisson))
## (Intercept) x
## 1.363235 1.650778
model_poisson$deviance / df.residual(model_poisson)
## [1] 1.083182
Jika nilai deviance/df > 1.5, kemungkinan terjadi overdispersion.
Visualisasi
data_poisson$pred <- predict(model_poisson, type = "response")
plot(data_poisson$x, data_poisson$y, pch = 16, col = "gray",
main = "Data dan Hasil Prediksi",
xlab = "x", ylab = "y")
ord <- order(data_poisson$x)
lines(data_poisson$x[ord], data_poisson$pred[ord], col = "blue", lwd = 2)
Ekspektasi
Jika diturunkan berdasarkan fungsi momen:
\[ \mathbb{E}(Y) = \int y f(y; \theta) dy = \mu \]
Untuk keluarga eksponensial:
\[ \log f_Y(y; \theta) = a(y) + b(\theta) y + c(\theta) \]
atau:
\[ \log f_Y(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) = \dot{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)
Implementasi Newton-Raphson
Statistik score ke-j:
\[ 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} = \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)}) \]
Diagnostik digunakan untuk mengevaluasi apakah model sudah tepat.
Statistik Devian
Devian adalah:
\[ 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 \frac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i} \]
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_i(\boldsymbol{x}) = \frac{\exp(\beta_0 + \beta_1 x_i)}{1 + \exp(\beta_0 + \beta_1 x_i)} \]
Log-likelihood untuk n observasi:
\[ \ell(\boldsymbol{\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(-\boldsymbol{x}_i^\top \boldsymbol{\beta})} \]
Log-likelihood untuk n observasi:
\[ \ell(\boldsymbol{\beta}) = \sum_{i=1}^n \left[ y_i \log(\pi_i) + (1 - y_i) \log(1 - \pi_i) \right] \]
Langkah-Langkah Newton-Raphson
\[ U(\beta) = \frac{\partial \ell(\beta)}{\partial \beta} = \mathbf{X}^\top (y - \pi) \]
\[ H(\beta) = -\mathbf{X}^\top \mathbf{W} \mathbf{X}, \quad \text{dengan } \mathbf{W} = \text{diag}(\pi_i(1 - \pi_i)) \]
\[ \beta^{(t+1)} = \beta^{(t)} + (\mathbf{X}^\top \mathbf{W}^{(t)} \mathbf{X})^{-1} \mathbf{X}^\top (y - \pi^{(t)}) \]
Estimasi MLE
# Data simulasi
set.seed(1)
x <- rnorm(100)
beta_true <- c(-1.5, 2.2)
X <- cbind(1, x)
eta <- X %*% beta_true
p <- 1 / (1 + exp(-eta))
y <- rbinom(100, 1, p)
# Iterasi Newton-Raphson
beta <- c(0, 0)
tol <- 1e-6
max_iter <- 150
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)) < 0.0001) break
beta <- beta_new
}
beta
## [,1]
## -2.295432
## x 2.561043
Uji Wald
# Data simulasi
set.seed(123)
x <- rnorm(150) # Sample size modified
beta_true <- c(-0.8, 1.5) # Parameters modified
X <- cbind(1, x)
eta <- X %*% beta_true
p <- 1 / (1 + exp(-eta))
y <- rbinom(150, 1, p)
# Fit model
model <- glm(y ~ x, family = binomial)
beta_hat <- coef(model)["x"]
se_beta <- summary(model)$coefficients["x", "Std. Error"] # Fixed space
Z <- beta_hat / se_beta
Wald_stat <- Z^2
p_value <- 1 - pchisq(Wald_stat, df = 1)
cat("Wald Statistic:", round(Wald_stat, 3),
"\np-value:", format.pval(p_value, digits = 3))
## Wald Statistic: 27.873
## p-value: 1.3e-07
Jika p-value < alpha (0.05), maka Terima Hipotesis nol.
Model regresi Poisson digunakan untuk memodelkan data count (jumlah kejadian) di mana variabel respon mengikuti distribusi Poisson. Estimasi dilakukan dengan Maximum Likelihood Estimation (MLE), dan inferensi dilakukan dengan uji Wald dan Likelihood Ratio Test.
Model Regresi Poisson
Fungsi distribusi Poisson:
\[ P(Y_i = y_i) = \frac{e^{-\lambda_i} \lambda_i^{y_i}}{y_i!} \]
Model regresi Poisson:
\[ \log(\lambda_i) = \mathbf{x}_i^\top \boldsymbol{\beta} \]
Estimasi Parameter (MLE)
Fungsi log-likelihood:
\[ \ell(\boldsymbol{\beta}) = \sum_{i=1}^n \left[ y_i \log(\lambda_i) - \lambda_i - \log(y_i!) \right] \]
Dengan:
\[ \lambda_i = \exp(\mathbf{x}_i^\top \boldsymbol{\beta}) \]
Estimasi menggunakan 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 optimisasi eksplisit.
Tahap 1: Definisikan model regresi Poisson:
\[ \log(\lambda_i) = \eta_i = \mathbf{x}_i^\top \boldsymbol{\beta} \]
sehingga:
\[ \lambda_i = \exp(\eta_i) \]
Tahap 2: Mencari log-likelihood yang dimaksimalkan:
\[ \ell(\boldsymbol{\beta}) = \sum_{i=1}^n \left[ y_i \log(\lambda_i) - \lambda_i - \log(y_i!) \right] \]
Tahap 3: Formulasi iteratif IRLS:
\[ \boldsymbol{\beta}^{(t+1)} = (\mathbf{X}^\top \mathbf{W} \mathbf{X})^{-1} \mathbf{X}^\top \mathbf{W} \mathbf{z} \]
Dengan:
dan
\[ \eta_i = \log(\lambda_i) = \mathbf{x}_i^\top \boldsymbol{\beta} \]
fEstimasi IRLS Manual
# Fungsi untuk estimasi manual regresi Poisson
poisson_irls <- function(X, y, max_iter = 50, tol = 1e-6) {
beta <- rep(0, ncol(X)) # Inisialisasi parameter
for (iter in 1:max_iter) {
eta <- X %*% beta
lambda <- exp(eta)
W <- diag(as.vector(lambda)) # Matriks bobot
z <- eta + (y - lambda) / lambda # Variabel respon adjusted
beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z
if (max(abs(beta_new - beta)) < tol) break
beta <- beta_new
}
list(coefficients = beta, iterations = iter)
}
##Contoh
set.seed(123)
n <- 200
X <- cbind(1, rnorm(n)) # Matriks desain dengan intercept
beta_true <- c(0.5, 0.8)
y <- rpois(n, exp(X %*% beta_true))
# Estimasi parameter
result <- poisson_irls(X, y)
cat("Estimasi parameter:\n")
## Estimasi parameter:
print(result$coefficients)
## [,1]
## [1,] 0.5095712
## [2,] 0.7923024
cat("Jumlah iterasi:", result$iterations)
## Jumlah iterasi: 9
Inferensi Statistik Uji Wald
# Menggunakan model dari glm() untuk contoh
model <- glm(y ~ X[,-1], family = poisson)
# Statistik Wald
wald_test <- function(model) {
summ <- summary(model)
coefs <- summ$coefficients
z <- coefs[, "Estimate"] / coefs[, "Std. Error"]
p_value <- 2 * pnorm(-abs(z))
data.frame(
Estimate = coefs[, "Estimate"],
Std.Error = coefs[, "Std. Error"],
z.value = z,
p.value = p_value
)
}
wald_results <- wald_test(model)
print(wald_results)
## Estimate Std.Error z.value p.value
## (Intercept) 0.5095711 0.05885608 8.657918 4.804596e-18
## X[, -1] 0.7923019 0.04433252 17.871800 1.955815e-71
Evaluasi Model
# AIC dan BIC
cat("AIC:", AIC(model), "\n")
## AIC: 635.3746
cat("BIC:", BIC(model), "\n")
## BIC: 641.9713
cat("Devians:", deviance(model))
## Devians: 216.0077
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: jeis pekerjaan (Freelance, Employee). Dalam analisis regresi
logistik, data nominal diubah menjadi variabel dummy (misal: Freelance =
1, selain itu = 0).
Ordinal: Variabel yang memiliki urutan logis
antar kategori, tetapi jarak antar kategori belum tentu sama.
Contoh: tingkat kepuasan pelanggan (Sangat Tidak Puas < Tidak Puas
< Netral < Puas < Sangat Puas). Dalam analisis, variabel
ordinal bisa diperlakukan:
Rasio: Variabel numerik kontinu yang memiliki
nol absolut dan rasio bermakna.
Contoh: jumlah jam bekerja per minggu. Data ini dapat langsung digunakan
dalam model tanpa transformasi khusus.
# Simulasi variabel
n <- 500
set.seed(123)
job_type <- sample(c("Freelance", "Employee"), n, replace = TRUE) # nominal
satisfaction <- sample(c("Very Unsatisfied", "Unsatisfied", "Neutral", "Satisfied", "Very Satisfied"),
n, replace = TRUE, prob = c(0.1, 0.2, 0.3, 0.25, 0.15)) # ordinal
work_hours <- rnorm(n, mean = 40, sd = 10) # rasio
# Logit model: asumsi dummy untuk Red, urutan satisfaction sebagai ordinal, dan jam kerja numerik
logit_p <- -1.5 +
0.6 * (job_type == "Freelance") +
0.7 * as.numeric(factor(satisfaction,
levels = c("Very Unsatisfied", "Unsatisfied", "Neutral", "Satisfied", "Very Satisfied"),
ordered = TRUE)) +
0.04 * work_hours
p <- 1 / (1 + exp(-logit_p))
set.seed(57)
success <- rbinom(n, 1, p)
library(tibble)
## Warning: package 'tibble' was built under R version 4.3.3
sim_data <- tibble(success, job_type, satisfaction, work_hours)
head(sim_data)
Interpretasi: Dataset berisikan 500 observasi dengan variabel warna mobil, kepuasan, jam kerja, dan status keberhasilan.
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.3.3
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
sim_data %>%
group_by(success) %>%
summarise(
Jumlah = n(),
Rata2_JamKerja = mean(work_hours),
.groups = "drop"
)
Interpretasi: Distribusi jumlah sukses dan tidak sukses adalah 445 dan 55. Rata-rata jam kerja pada masing-masing grup adalah 35.28656 dan 40.55633.
sim_data_nominal <- sim_data %>%
mutate(
satisfaction = factor(satisfaction, levels = c("Very Unsatisfied", "Unsatisfied", "Neutral", "Satisfied", "Very Satisfied"))
)
model_nominal <- glm(success ~ job_type + satisfaction + work_hours, data = sim_data_nominal, family = binomial)
summary(model_nominal)
##
## Call:
## glm(formula = success ~ job_type + satisfaction + work_hours,
## family = binomial, data = sim_data_nominal)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.78678 0.73609 -2.427 0.015207 *
## job_typeFreelance 0.83875 0.31746 2.642 0.008240 **
## satisfactionUnsatisfied 0.56447 0.45251 1.247 0.212240
## satisfactionNeutral 1.25094 0.46212 2.707 0.006790 **
## satisfactionSatisfied 1.90797 0.52677 3.622 0.000292 ***
## satisfactionVery Satisfied 3.63136 1.07618 3.374 0.000740 ***
## work_hours 0.06026 0.01584 3.804 0.000143 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 346.52 on 499 degrees of freedom
## Residual deviance: 294.61 on 493 degrees of freedom
## AIC: 308.61
##
## Number of Fisher Scoring iterations: 7
Keterangan: Model ini menggunakan satisfaction sebagai variabel dummy, baseline ““Very Unsatisfied”. Setiap koefisien membandingkan kategori terhadap baseline.
Interpretasi Koefisisen
Intercept (-1.78678)
car_colorBlue (0.83875)
satisfactionUnsatisfied (0.56447)
satisfactionNeutral (1.25094)
satisfactionSatisfied (1.90797)
work_hours (0.06026)
Interpretasi Goodness-of-Fit
Null Deviance (346.52): Deviance model tanpa prediktor
Residual deviance (294.61): Deviance model dengan prediktor
AIC (308.61): Digunakan untuk membandingkan model, semakin kecil AIC maka semakin baik model dalam menyeimbangkan goodness-of-fit dan kompleksitas.
Signifikansi Model
Variabel satisfaction secara keseluruhan signifikan dalam meningkatkan peluang sukses. Terutama kategori Neutral, Satisfied, dan Very Satisfied, yang memiliki pengaruh signifikan (p < 0.01).
Warna mobil Blue juga menunjukkan pengaruh signifikan terhadap peluang sukses (p = 0.0082), dengan peluang sukses sekitar 2.3 kali lebih tinggi dibanding warna referensi.
Jumlah jam kerja (work_hours) signifikan secara statistik (p < 0.001), artinya semakin banyak jam kerja, semakin tinggi peluang sukses individu.
Kategori satisfaction Unsatisfied memiliki efek positif namun tidak signifikan (p = 0.2122), sehingga tidak cukup bukti bahwa kategori ini berbeda secara signifikan dari Very Unsatisfied.
Kesimpulan Praktis
Tingkat kepuasan individu merupakan prediktor paling kuat terhadap keberhasilan. Semakin tinggi tingkat kepuasan, semakin besar peluang sukses, dengan peningkatan hingga 38 kali lipat untuk kategori Very Satisfied dibanding Very Unsatisfied.
Warna mobil biru secara mengejutkan memiliki hubungan positif terhadap keberhasilan, namun ini bisa jadi artefak dari data simulasi. Interpretasi ini perlu dikaji lebih lanjut dengan variabel pengganti yang lebih logis.
Jam kerja memberikan kontribusi positif yang signifikan terhadap keberhasilan, dengan tambahan 1 jam kerja meningkatkan peluang sukses sebesar 6.2%.
Meskipun model ini sudah menjelaskan variabilitas data lebih baik daripada model null (lihat deviance dan AIC), masih ada ruang untuk perbaikan, misalnya dengan menambahkan variabel yang lebih bermakna seperti jenis pekerjaan atau motivasi pribadi.
library(dplyr)
sim_data_numeric <- sim_data %>%
mutate(
satisfaction_numeric = case_when(
satisfaction == "Very Unsatisfied" ~ 1,
satisfaction == "Unsatisfied" ~ 2,
satisfaction == "Neutral" ~ 3,
satisfaction == "Satisfied" ~ 4,
satisfaction == "Very Satisfied" ~ 5
)
)
model_numeric <- glm(success ~ job_type + satisfaction_numeric + work_hours, data = sim_data_numeric, family = binomial)
summary(model_numeric)
##
## Call:
## glm(formula = success ~ job_type + satisfaction_numeric + work_hours,
## family = binomial, data = sim_data_numeric)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.68665 0.75758 -3.546 0.000391 ***
## job_typeFreelance 0.86532 0.31581 2.740 0.006144 **
## satisfaction_numeric 0.72901 0.14034 5.195 2.05e-07 ***
## work_hours 0.06099 0.01590 3.836 0.000125 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 346.52 on 499 degrees of freedom
## Residual deviance: 295.91 on 496 degrees of freedom
## AIC: 303.91
##
## Number of Fisher Scoring iterations: 6
Interpretasi: Model menggunakan satisfaction sebagai angka bertingkat 1-5. Koefisien mengukur kenaikan satu tingkat kepuasan terhadap peluang sukses.
Perbandingan Model
list(
AIC_Nominal = AIC(model_nominal),
AIC_Numeric = AIC(model_numeric)
)
## $AIC_Nominal
## [1] 308.6077
##
## $AIC_Numeric
## [1] 303.9145
Keterangan: Model dengan nilai AIC yang lebih kecil cenderung lebih disukai.
Interpretasi Koefisien Intercept (-2.68665) Ini adalah log-odds dasar untuk individu dengan tipe pekerjaan Employee, kepuasan kerja terendah (Very Unsatisfied = 1), dan jam kerja = 0. Tidak dapat diinterpretasikan secara praktis karena tidak realistis (jam kerja = 0), namun tetap menjadi baseline referensi.
job_typeFreelance (0.86532) Individu dengan tipe pekerjaan Freelance memiliki log-odds peluang sukses yang lebih tinggi sebesar 0.865 dibanding Employee.
p-value = 0.006144 (signifikan)
Odds Ratio = exp(0.86532) = 2.376 → Peluang sukses freelancer sekitar 2.4 kali lebih tinggi dibandingkan pegawai.
satisfaction_numeric (0.72901) Setiap peningkatan satu tingkat kepuasan kerja (misalnya dari Unsatisfied ke Neutral) meningkatkan log-odds sukses sebesar 0.729.
p-value < 0.001 (sangat signifikan)
Odds Ratio = exp(0.72901) = 2.073 → Setiap kenaikan tingkat kepuasan kerja meningkatkan peluang sukses sebesar 107%.
work_hours (0.06099) Setiap penambahan satu jam kerja per minggu meningkatkan log-odds sukses sebesar 0.061.
p-value = 0.000125 (signifikan)
Odds Ratio = exp(0.06099) = 1.063 → Setiap tambahan satu jam kerja per minggu meningkatkan peluang sukses sekitar 6.3%.
Interpretasi Goodness-of-Fit - Null deviance (346.52): Deviance model tanpa prediktor.
Residual deviance (295.91): Deviance model dengan prediktor.→ Penurunan deviance menunjukkan bahwa prediktor memberikan tambahan informasi yang signifikan.
AIC (303.91): Digunakan untuk membandingkan model — semakin kecil, semakin baik.
Signifikansi Model Variabel satisfaction_numeric, work_hours, dan job_typeFreelance semuanya signifikan dalam meningkatkan peluang sukses.
Model secara umum baik dalam menjelaskan variabilitas respons (success).
Kesimpulan - Kepuasan kerja dan jam kerja memiliki pengaruh positif terhadap peluang sukses.
Pekerja freelance memiliki peluang lebih tinggi dibandingkan pekerja tetap.
Model yang dibangun secara statistik signifikan dan mampu menjelaskan variansi data.
Goodeness-of-Fit Model
nullmod <- glm(success ~ 1, data = sim_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.1497989 (df=7)
##
## $McFadden_R2_Numeric
## 'log Lik.' 0.1460278 (df=4)
Interpretasi: McFadden’s R^2 mengukur goodness-of-fit. Semakin besar nilainya, semakin baik model memprediksi dibandingkan model null.
Visualisasi Prediksi
library(ggplot2)
sim_data_nominal <- sim_data_nominal %>% mutate(predicted = predict(model_nominal, type = "response"))
sim_data_numeric <- sim_data_numeric %>% mutate(predicted = predict(model_numeric, type = "response"))
# Plot untuk model nominal
sim_data_nominal %>%
ggplot(aes(x = work_hours, y = predicted, color = satisfaction)) +
geom_point(alpha = 0.6) +
labs(title = "Prediksi Probabilitas (Ordinal sebagai Nominal)", x = "work hours", y = "Prediksi Probabiliti")
theme_minimal()
## List of 136
## $ line :List of 6
## ..$ colour : chr "black"
## ..$ linewidth : num 0.5
## ..$ linetype : num 1
## ..$ lineend : chr "butt"
## ..$ arrow : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_line" "element"
## $ rect :List of 5
## ..$ fill : chr "white"
## ..$ colour : chr "black"
## ..$ linewidth : num 0.5
## ..$ linetype : num 1
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ text :List of 11
## ..$ family : chr ""
## ..$ face : chr "plain"
## ..$ colour : chr "black"
## ..$ size : num 11
## ..$ hjust : num 0.5
## ..$ vjust : num 0.5
## ..$ angle : num 0
## ..$ lineheight : num 0.9
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ title : NULL
## $ aspect.ratio : NULL
## $ axis.title : NULL
## $ axis.title.x :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 2.75points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.x.top :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 0
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 2.75points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.x.bottom : NULL
## $ axis.title.y :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : num 90
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 2.75points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.y.left : NULL
## $ axis.title.y.right :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : num -90
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 2.75points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : chr "grey30"
## ..$ size : 'rel' num 0.8
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 2.2points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x.top :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 0
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 2.2points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x.bottom : NULL
## $ axis.text.y :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 1
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 2.2points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.y.left : NULL
## $ axis.text.y.right :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 2.2points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.theta : NULL
## $ axis.text.r :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0.5
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 2.2points 0points 2.2points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.ticks : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ axis.ticks.x : NULL
## $ axis.ticks.x.top : NULL
## $ axis.ticks.x.bottom : NULL
## $ axis.ticks.y : NULL
## $ axis.ticks.y.left : NULL
## $ axis.ticks.y.right : NULL
## $ axis.ticks.theta : NULL
## $ axis.ticks.r : NULL
## $ axis.minor.ticks.x.top : NULL
## $ axis.minor.ticks.x.bottom : NULL
## $ axis.minor.ticks.y.left : NULL
## $ axis.minor.ticks.y.right : NULL
## $ axis.minor.ticks.theta : NULL
## $ axis.minor.ticks.r : NULL
## $ axis.ticks.length : 'simpleUnit' num 2.75points
## ..- attr(*, "unit")= int 8
## $ axis.ticks.length.x : NULL
## $ axis.ticks.length.x.top : NULL
## $ axis.ticks.length.x.bottom : NULL
## $ axis.ticks.length.y : NULL
## $ axis.ticks.length.y.left : NULL
## $ axis.ticks.length.y.right : NULL
## $ axis.ticks.length.theta : NULL
## $ axis.ticks.length.r : NULL
## $ axis.minor.ticks.length : 'rel' num 0.75
## $ axis.minor.ticks.length.x : NULL
## $ axis.minor.ticks.length.x.top : NULL
## $ axis.minor.ticks.length.x.bottom: NULL
## $ axis.minor.ticks.length.y : NULL
## $ axis.minor.ticks.length.y.left : NULL
## $ axis.minor.ticks.length.y.right : NULL
## $ axis.minor.ticks.length.theta : NULL
## $ axis.minor.ticks.length.r : NULL
## $ axis.line : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ axis.line.x : NULL
## $ axis.line.x.top : NULL
## $ axis.line.x.bottom : NULL
## $ axis.line.y : NULL
## $ axis.line.y.left : NULL
## $ axis.line.y.right : NULL
## $ axis.line.theta : NULL
## $ axis.line.r : NULL
## $ legend.background : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ legend.margin : 'margin' num [1:4] 5.5points 5.5points 5.5points 5.5points
## ..- attr(*, "unit")= int 8
## $ legend.spacing : 'simpleUnit' num 11points
## ..- attr(*, "unit")= int 8
## $ legend.spacing.x : NULL
## $ legend.spacing.y : NULL
## $ legend.key : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ legend.key.size : 'simpleUnit' num 1.2lines
## ..- attr(*, "unit")= int 3
## $ legend.key.height : NULL
## $ legend.key.width : NULL
## $ legend.key.spacing : 'simpleUnit' num 5.5points
## ..- attr(*, "unit")= int 8
## $ legend.key.spacing.x : NULL
## $ legend.key.spacing.y : NULL
## $ legend.frame : NULL
## $ legend.ticks : NULL
## $ legend.ticks.length : 'rel' num 0.2
## $ legend.axis.line : NULL
## $ legend.text :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : 'rel' num 0.8
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ legend.text.position : NULL
## $ legend.title :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ legend.title.position : NULL
## $ legend.position : chr "right"
## $ legend.position.inside : NULL
## $ legend.direction : NULL
## $ legend.byrow : NULL
## $ legend.justification : chr "center"
## $ legend.justification.top : NULL
## $ legend.justification.bottom : NULL
## $ legend.justification.left : NULL
## $ legend.justification.right : NULL
## $ legend.justification.inside : NULL
## $ legend.location : NULL
## $ legend.box : NULL
## $ legend.box.just : NULL
## $ legend.box.margin : 'margin' num [1:4] 0cm 0cm 0cm 0cm
## ..- attr(*, "unit")= int 1
## $ legend.box.background : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ legend.box.spacing : 'simpleUnit' num 11points
## ..- attr(*, "unit")= int 8
## [list output truncated]
## - attr(*, "class")= chr [1:2] "theme" "gg"
## - attr(*, "complete")= logi TRUE
## - attr(*, "validate")= logi TRUE
# Plot untuk model numeric
sim_data_numeric %>%
ggplot(aes(x = work_hours, y = predicted, color = as.factor(satisfaction_numeric))) +
geom_point(alpha = 0.6) +
labs(title = "Prediksi Probabilitas (Ordinal sebagai Numeric)", x = "Work Hours", y = "Prediksi Probabiliti")
theme_minimal()
## List of 136
## $ line :List of 6
## ..$ colour : chr "black"
## ..$ linewidth : num 0.5
## ..$ linetype : num 1
## ..$ lineend : chr "butt"
## ..$ arrow : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_line" "element"
## $ rect :List of 5
## ..$ fill : chr "white"
## ..$ colour : chr "black"
## ..$ linewidth : num 0.5
## ..$ linetype : num 1
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_rect" "element"
## $ text :List of 11
## ..$ family : chr ""
## ..$ face : chr "plain"
## ..$ colour : chr "black"
## ..$ size : num 11
## ..$ hjust : num 0.5
## ..$ vjust : num 0.5
## ..$ angle : num 0
## ..$ lineheight : num 0.9
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : logi FALSE
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ title : NULL
## $ aspect.ratio : NULL
## $ axis.title : NULL
## $ axis.title.x :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 2.75points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.x.top :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 0
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 2.75points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.x.bottom : NULL
## $ axis.title.y :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : num 90
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 2.75points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.title.y.left : NULL
## $ axis.title.y.right :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : num -90
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 2.75points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : chr "grey30"
## ..$ size : 'rel' num 0.8
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 1
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 2.2points 0points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x.top :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : NULL
## ..$ vjust : num 0
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 2.2points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.x.bottom : NULL
## $ axis.text.y :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 1
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 2.2points 0points 0points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.y.left : NULL
## $ axis.text.y.right :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 0points 0points 2.2points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.text.theta : NULL
## $ axis.text.r :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0.5
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : 'margin' num [1:4] 0points 2.2points 0points 2.2points
## .. ..- attr(*, "unit")= int 8
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ axis.ticks : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ axis.ticks.x : NULL
## $ axis.ticks.x.top : NULL
## $ axis.ticks.x.bottom : NULL
## $ axis.ticks.y : NULL
## $ axis.ticks.y.left : NULL
## $ axis.ticks.y.right : NULL
## $ axis.ticks.theta : NULL
## $ axis.ticks.r : NULL
## $ axis.minor.ticks.x.top : NULL
## $ axis.minor.ticks.x.bottom : NULL
## $ axis.minor.ticks.y.left : NULL
## $ axis.minor.ticks.y.right : NULL
## $ axis.minor.ticks.theta : NULL
## $ axis.minor.ticks.r : NULL
## $ axis.ticks.length : 'simpleUnit' num 2.75points
## ..- attr(*, "unit")= int 8
## $ axis.ticks.length.x : NULL
## $ axis.ticks.length.x.top : NULL
## $ axis.ticks.length.x.bottom : NULL
## $ axis.ticks.length.y : NULL
## $ axis.ticks.length.y.left : NULL
## $ axis.ticks.length.y.right : NULL
## $ axis.ticks.length.theta : NULL
## $ axis.ticks.length.r : NULL
## $ axis.minor.ticks.length : 'rel' num 0.75
## $ axis.minor.ticks.length.x : NULL
## $ axis.minor.ticks.length.x.top : NULL
## $ axis.minor.ticks.length.x.bottom: NULL
## $ axis.minor.ticks.length.y : NULL
## $ axis.minor.ticks.length.y.left : NULL
## $ axis.minor.ticks.length.y.right : NULL
## $ axis.minor.ticks.length.theta : NULL
## $ axis.minor.ticks.length.r : NULL
## $ axis.line : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ axis.line.x : NULL
## $ axis.line.x.top : NULL
## $ axis.line.x.bottom : NULL
## $ axis.line.y : NULL
## $ axis.line.y.left : NULL
## $ axis.line.y.right : NULL
## $ axis.line.theta : NULL
## $ axis.line.r : NULL
## $ legend.background : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ legend.margin : 'margin' num [1:4] 5.5points 5.5points 5.5points 5.5points
## ..- attr(*, "unit")= int 8
## $ legend.spacing : 'simpleUnit' num 11points
## ..- attr(*, "unit")= int 8
## $ legend.spacing.x : NULL
## $ legend.spacing.y : NULL
## $ legend.key : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ legend.key.size : 'simpleUnit' num 1.2lines
## ..- attr(*, "unit")= int 3
## $ legend.key.height : NULL
## $ legend.key.width : NULL
## $ legend.key.spacing : 'simpleUnit' num 5.5points
## ..- attr(*, "unit")= int 8
## $ legend.key.spacing.x : NULL
## $ legend.key.spacing.y : NULL
## $ legend.frame : NULL
## $ legend.ticks : NULL
## $ legend.ticks.length : 'rel' num 0.2
## $ legend.axis.line : NULL
## $ legend.text :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : 'rel' num 0.8
## ..$ hjust : NULL
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ legend.text.position : NULL
## $ legend.title :List of 11
## ..$ family : NULL
## ..$ face : NULL
## ..$ colour : NULL
## ..$ size : NULL
## ..$ hjust : num 0
## ..$ vjust : NULL
## ..$ angle : NULL
## ..$ lineheight : NULL
## ..$ margin : NULL
## ..$ debug : NULL
## ..$ inherit.blank: logi TRUE
## ..- attr(*, "class")= chr [1:2] "element_text" "element"
## $ legend.title.position : NULL
## $ legend.position : chr "right"
## $ legend.position.inside : NULL
## $ legend.direction : NULL
## $ legend.byrow : NULL
## $ legend.justification : chr "center"
## $ legend.justification.top : NULL
## $ legend.justification.bottom : NULL
## $ legend.justification.left : NULL
## $ legend.justification.right : NULL
## $ legend.justification.inside : NULL
## $ legend.location : NULL
## $ legend.box : NULL
## $ legend.box.just : NULL
## $ legend.box.margin : 'margin' num [1:4] 0cm 0cm 0cm 0cm
## ..- attr(*, "unit")= int 1
## $ legend.box.background : list()
## ..- attr(*, "class")= chr [1:2] "element_blank" "element"
## $ legend.box.spacing : 'simpleUnit' num 11points
## ..- attr(*, "unit")= int 8
## [list output truncated]
## - attr(*, "class")= chr [1:2] "theme" "gg"
## - attr(*, "complete")= logi TRUE
## - attr(*, "validate")= logi TRUE
Interpretasi: Visualisasi hubungan antara pendapatan dan probabilitas sukses berdasarkan tingkat kepuasan dengan dua pendekatan ordinal. Secara numerik, memperhatikan urutan dari Very Unsatisfied hingga Very Satisfied. Namun, secara nominal tidak memperhatikan urutan tingkat kepuasan nya.
Ringkasan Koefisien Model
library(broom)
## Warning: package 'broom' was built under R version 4.3.3
library(knitr)
## Warning: package 'knitr' was built under R version 4.3.3
library(kableExtra)
## Warning: package 'kableExtra' was built under R version 4.3.3
##
## Attaching package: 'kableExtra'
## The following object is masked from 'package:dplyr':
##
## group_rows
# Ringkasan model nominal
tidy(model_nominal) %>%
kable(format = "latex", booktabs = TRUE,
caption = "Ringkasan Koefisien Model dengan Ordinal sebagai Nominal") %>%
kable_styling(latex_options = c("hold_position","striped"))
# Ringkasan model numeric
tidy(model_numeric) %>%
kable(format = "latex", booktabs = TRUE, caption = "Ringkasan Koefisien Model dengan Ordinal sebagai Numeric") %>% kable_styling(latex_options = c("hold_position", "striped"))
Interpretasi: Tabel memperlihatkan estimasi koefisien, standar error, nilai z, dan p-value. Koefisien positif meningkatkan log-odds keberhasilan, koefisien negatif akan menurunkan log-odds keberhasilan.
Kesimpulan
Job Type: Freelance memiliki peluang sukses lebih tinggi dibandingkan dengan employee
Satisfaction:
Jika dianggap dummy (nominal) : tiap tingkat akan dibandingkan langsung dengan tingkat Very Unsatisfied
Jika dianggap numeric : tiap kenaikan satu tingkat kepuasan meningkatkan peluang sukses.
Work Hours: Waktu bekerja yang lebih tinggi meningkatkan peluang sukses.
Pada akhirnya, model dengan perlakuan ordinal sebagai numerik akan memiliki interpretasi lebih sederhana jika asumsi jarak antar tingkatannya konsisten.
Dalam analisis regresi logistik, pemilihan model yang tepat sangat penting untuk memperoleh model yang dapat memprediksi probabilitas kejadian suatu peristiwa (respon biner) secara akurat. Dua pendekatan utama yang sering digunakan dalam membangun model adalah pendekatan Konfirmatori dan Eksploratori.
1. Konfirmatori (Pendekatan Konfirmatori)
Pendekatan ini digunakan apabila peneliti sudah memiliki teori atau hipotesis yang jelas mengenai hubungan antara variabel prediktor dengan variabel respon.
Ciri-ciri:
Contoh penggunaan:
Misalnya, berdasarkan teori diasumsikan bahwa faktor x1 dan x2 mempengaruhi kemungkinan seseorang membeli produk. Maka, model regresi logistik akan langsung mengikutkan x1 dan x2 untuk diuji apakah kontribusinya signifikan.
2. Eksploratori (Pendekatan Eksploratori)
Pendekatan ini digunakan ketika peneliti belum memiliki teori yang jelas, dan bertujuan untuk mengeksplorasi hubungan antara variabel.
Ciri-ciri:
Proses pemilihan model:
Tujuan:
Mendapatkan model yang sederhana namun cukup akurat dalam memprediksi
nilai respon.
Pemilihan antara pendekatan Konfirmatori atau Eksploratori bergantung
pada tujuan penelitian. Jika terdapat teori yang jelas dan mendasar,
maka pendekatan konfirmatori lebih sesuai. Namun, jika belum ada teori
yang kuat, pendekatan eksploratori bisa membantu menemukan hubungan yang
relevan untuk diuji lebih lanjut dalam penelitian selanjutnya.
Simulasi Data
# Membangun Model Regresi Logistik
library(knitr)
library(dplyr)
library(ggplot2)
library(MASS)
##
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
##
## select
library(caret)
## Warning: package 'caret' was built under R version 4.3.3
## Loading required package: lattice
## Warning: package 'lattice' was built under R version 4.3.3
##
## Attaching package: 'caret'
## The following objects are masked from 'package:DescTools':
##
## MAE, RMSE
library(pROC)
## Warning: package 'pROC' was built under R version 4.3.3
## Type 'citation("pROC")' for a citation.
##
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
##
## cov, smooth, var
library(DescTools)
## Simulasi Data
set.seed(57)
n <- 200
x1 <- rnorm(n)
x2 <- rbinom(n, 1, 0.5)
x3 <- rnorm(n)
lin_pred <- -0.5 + 1.5 * x1 + 0.4 * x2 - 0.9 * x3
p <- 1 / (1 + exp(-lin_pred))
y <- rbinom(n, 1, p)
df <- data.frame(y = as.factor(y), x1, x2, x3)
head(df)
Pemilihan Model Model Full
### 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.5254 0.2708 -1.941 0.0523 .
## x1 1.4894 0.2315 6.433 1.25e-10 ***
## x2 0.4708 0.3692 1.275 0.2023
## x3 -0.8609 0.1982 -4.343 1.41e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 276.28 on 199 degrees of freedom
## Residual deviance: 190.58 on 196 degrees of freedom
## AIC: 198.58
##
## Number of Fisher Scoring iterations: 5
### Metode Forward, Backward, dan Stepwise
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)
Interpretasi: Model dengan AIC lebih rendah umumnya lebih disukai untuk dipilih menjadi model terbaik.
## Evaluasi Model : ROC dan AUC
pred_prob <- predict(step_both, type = "response")
roc_obj <- roc(df$y, pred_prob)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
plot(roc_obj, main = "Kurva ROC", col = "red")
auc(roc_obj)
## Area under the curve: 0.8444
## Pseudo R-Squared
PseudoR2(step_both, which = c("CoxSnell", "Nagelkerke", "McFadden"))
## CoxSnell Nagelkerke McFadden
## 0.3431234 0.4582487 0.3042290
Interpretasi: 1. Cox & Snell (0.2386): Nilai ini menunjukkan bahwa model menjelaskan sekitar 23.86% variasi dalam data. Namun, nilai ini memiliki batas atas yang lebih rendah dari 1, sehingga interpretasinya terbatas.
Nagelkerke (0.3294): Merupakan penyesuaian dari Cox & Snell agar skalanya berkisar antara 0 hingga 1. Dengan nilai 32.94%, ini menunjukkan bahwa model menjelaskan proporsi variasi yang cukup baik dan lebih mudah diinterpretasikan dibandingkan Cox & Snell.
McFadden (0.2115): Mengindikasikan bahwa model memiliki kemampuan prediktif yang cukup layak. Nilai antara 0.2 hingga 0.4 umumnya dianggap sebagai indikasi model yang baik dalam konteks regresi logistik.
## Tabel Klasifikasi dan Evaluasi
pred_class <- ifelse(pred_prob >= 0.5, 1, 0)
conf_matrix <- confusionMatrix(factor(pred_class), df$y, positive = "1")
conf_matrix
## Confusion Matrix and Statistics
##
## Reference
## Prediction 0 1
## 0 86 25
## 1 21 68
##
## Accuracy : 0.77
## 95% CI : (0.7054, 0.8264)
## No Information Rate : 0.535
## P-Value [Acc > NIR] : 5.176e-12
##
## Kappa : 0.5364
##
## Mcnemar's Test P-Value : 0.6583
##
## Sensitivity : 0.7312
## Specificity : 0.8037
## Pos Pred Value : 0.7640
## Neg Pred Value : 0.7748
## Prevalence : 0.4650
## Detection Rate : 0.3400
## Detection Prevalence : 0.4450
## Balanced Accuracy : 0.7675
##
## 'Positive' Class : 1
##
conf_matrix$byClass[c("Sensitivity", "Specificity")]
## Sensitivity Specificity
## 0.7311828 0.8037383
Model terbaik berdasarkan AIC adalah model step both dengan AIC paling minimum. Dengan nilai ROC menunjukan 0.8444 yang terbilang cukup baik. Nilai Pseudo R-Squared dan hasil klasifikasi model juga menunjukkan hasil yang cuku baik. Dimana klasifikasi memiliki akurasi sebesar 77%.
Tujuan Dalam sub-bab ini akan disajikan cara membandingkan model regresi logistik menggunakan ukuran Deviance, AIC, dan Likelihood-Ratio serta menjelaskan prinsip Parsimony dalam pemilihan model.
Simulasi Data
# Metode Perbandingan Model dalam Regresi Logistik
library(MASS)
library(broom)
library(DescTools)
## Simulasi Data
set.seed(57)
n <- 300
x1 <- rnorm(n)
x2 <- rbinom(n, 1, 0.5)
x3 <- rnorm(n)
lin_pred <- -1 + 1.5 * x1 - 0.7 * x2 + 0.9 * x3
p <- 1 / (1 + exp(-lin_pred))
y <- rbinom(n, 1, p)
data <- data.frame(y = as.factor(y), x1, x2, x3)
Pembuatan Model
## Pembuatan Model
model1 <- glm(y ~ x1, data = data, family = binomial)
model2 <- glm(y ~ x1 + x2, data = data, family = binomial)
model3 <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
Perbandingan AIC dan Deviance
## 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
Dari perbandingan nilai AIC dan Deviance dari 3 model, model terbaik adalah model dengan AIC dan Deviance paling kecil. Hal itu terdapat pada model dengan variabel paling banyak.
anova(model1, model2, test = "LRT")
Uji Likelihood-Rasio Test menunjukkan p-value 0.001768. P-value < 0.05 maka secara signifikan penambahan variabel x2 meningkatkan kecocokan pada model. x2 layak dipertahankan.
anova(model2, model3, test = "LRT")
Uji Likelihood-Rasio Test menunjukkan p-value 4.849e-11. P-value < 0.05 maka secara signifikan penambahan variabel x2 meningkatkan kecocokan pada model. x3 layak dipertahankan.
Dalam pemodelan, model yang kompleks umumnya menghasilkan nilai AIC dan deviance yang lebih kecil. Namun, model yang lebih sederhana lebih mudah untuk dipahami dan dijelaskan. Jika perbedaan AIC antar model tidak signifikan, maka model yang lebih sederhana sebaiknya dipilih. Model kompleks yang terlalu overfit bisa jadi tidak generalis.
Rumus dan Penjabaran
Rumus AIC
\[ \text{AIC} = -2(\log L - k) = -2 \log L + 2k \]
Penjelasan:
AIC (Akaike Information Criterion) digunakan untuk menilai seberapa baik
model sesuai dengan data (melalui likelihood) dengan penalti untuk
kompleksitas model (jumlah parameter, \(k\)).
Model dengan nilai AIC lebih kecil dianggap lebih baik.
Namun, jika perbedaannya tidak signifikan, pilihlah model yang
lebih sederhana, karena model yang terlalu kompleks rentan
terhadap overfitting.
Rumus Deviance
\[ D = -2 \log(L_{\text{model}}) - \log(L_{\text{model\_saturasi}}) \]
Penjelasan: Deviance adalah ukuran seberapa jauh model saat ini dari model yang sangat sesuai (model saturasi). Semakin kecil nilai deviance, semakin baik model tersebut dalam menjelaskan variasi dalam data.
Rumus Likelihood Ratio Test
\[ G^2 = -2 \log L_0 + 2 \log L_1 \]
Penjelasan: Statistik Likelihood Ratio Test
(LRT) digunakan untuk membandingkan dua model yang bersarang
(nested). Model yang lebih kompleks diuji apakah memberikan peningkatan
signifikan terhadap kecocokan model dibanding model yang lebih
sederhana.
Jika nilai \(G^2\) besar dan
p-value signifikan, maka model kompleks dianggap lebih
baik dari model sederhana.
## Evaluasi Tabel Klasifikasi dan Akurasi Model
pred_prob <- predict(model3, type = "response")
pred_class <- factor(ifelse(pred_prob >= 0.5, 1, 0))
conf_matrix <- confusionMatrix(pred_class, data$y, positive = "1")
conf_matrix
## Confusion Matrix and Statistics
##
## Reference
## Prediction 0 1
## 0 194 38
## 1 18 50
##
## Accuracy : 0.8133
## 95% CI : (0.7646, 0.8558)
## No Information Rate : 0.7067
## P-Value [Acc > NIR] : 1.589e-05
##
## Kappa : 0.5177
##
## Mcnemar's Test P-Value : 0.01112
##
## Sensitivity : 0.5682
## Specificity : 0.9151
## Pos Pred Value : 0.7353
## Neg Pred Value : 0.8362
## Prevalence : 0.2933
## Detection Rate : 0.1667
## Detection Prevalence : 0.2267
## Balanced Accuracy : 0.7416
##
## 'Positive' Class : 1
##
conf_matrix$byClass[c("Sensitivity", "Specificity")]
## Sensitivity Specificity
## 0.5681818 0.9150943
Kesimpulan - Deviance yang kecil menunjukkan kecocokan model yang lebih baik - AIC yang lebih rendah menunjukkan keseimbangan antara kecocokan dan kompleksitas - Likelihood-Ratio Test mengevaluasi apakah model kompleks secara signifikan lebih baik - Tabel Klasifikasi membantu menilai kinerja prediksi aktual dengan prediksi model - Prinsip Parsimony mengutamakan model sederhana jika performanya cenderung mirip.
Kurva ROC (Receiver Operating Characteristic) merupakan alat visual yang sering digunakan untuk mengevaluasi performa dari model klasifikasi biner. Kurva ini memperlihatkan hubungan antara tingkat True Positive Rate (Sensitivity) dan False Positive Rate (1 - Specificity) pada berbagai nilai ambang klasifikasi.
Perubahan nilai ambang klasifikasi (cut-off) dapat mempengaruhi banyaknya data yang diklasifikasikan sebagai positif:
Kurva ROC yang ideal memiliki ciri-ciri berikut:
Interpretasi dari nilai AUC adalah sebagai berikut:
Nilai AUC juga bisa diinterpretasikan sebagai probabilitas bahwa model dapat memberikan skor yang lebih tinggi kepada kasus positif dibandingkan kasus negatif.
Kurva ROC digunakan untuk:
Di dalam R, kurva ROC dapat dibuat menggunakan package
pROC.
library(pROC)
set.seed(57)
x1 <- rnorm(200)
x2 <- rbinom(200, 1, 0.5)
x3 <- rnorm(200)
lin_pred <- -1 + 1.3 * x1 - 0.5 * x2 + 0.7 * x3
p <- 1 / (1 + exp(-lin_pred))
y <- rbinom(200, 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.8389
thresholds <- seq(0.1, 0.9, by = 0.05)
results <- data.frame(Threshold = thresholds)
results$Sensitivity <- sapply(thresholds, function(t) {
pred_class <- ifelse(pred >= t, 1, 0)
cm <- table(Pred = pred_class, Obs = data$y)
TP <- cm["1", "1"]
FN <- cm["0", "1"]
TP / (TP + FN)
})
results$Specificity <- sapply(thresholds, function(t) {
pred_class <- ifelse(pred >= t, 1, 0)
cm <- table(Pred = pred_class, Obs = data$y)
TN <- cm["0", "0"]
FP <- cm["1", "0"]
TN / (TN + FP)
})
print(results)
## Threshold Sensitivity Specificity
## 1 0.10 0.9866667 0.336
## 2 0.15 0.9200000 0.448
## 3 0.20 0.8800000 0.504
## 4 0.25 0.8133333 0.600
## 5 0.30 0.8133333 0.672
## 6 0.35 0.8000000 0.704
## 7 0.40 0.7333333 0.736
## 8 0.45 0.7200000 0.824
## 9 0.50 0.6800000 0.864
## 10 0.55 0.6133333 0.888
## 11 0.60 0.5733333 0.904
## 12 0.65 0.5066667 0.928
## 13 0.70 0.4000000 0.984
## 14 0.75 0.3200000 0.992
## 15 0.80 0.2666667 0.992
## 16 0.85 0.1866667 1.000
## 17 0.90 0.0800000 1.000
Cut-off optimal dapat dipilih berdasarkan: - Nilai maksimum dari Sensitivity + Spesifisitas - Mempertimbangkan trade-off sesuai tujuan aplikasi
Penggunaan ROC cocok digunakan saat proporsi kelas seimbang. Namun, jika tidak precision-recall curve bisa lebih informatif.
Penjelasan Precision-Recall Curve Kurva Precision-Recall (PR) adalah alat evaluasi performa model klasifkasi, khususnya sangat berguna saat bekerja dengan data yang tidak seimbang (class imbalance).
Presisi (Precision): Merupakan rasio prediksi positif yang benar, yaitu berapa banyak dari seluruh prediksi positif yang benar-benar relevan. \[ \text{Precision} = \frac{TP}{TP + FP} \]
Recall (Sensitivitas): Merupakan rasio kasus positif yang berhasil diidentifikasi dengan benar sebagai positif. \[ \text{Recall} = \frac{TP}{TP + FN} \]
| Aspek | Kurva ROC | Kurva Precision-Recall |
|---|---|---|
| Fokus | Semua kelas | Kelas positif saja |
| Kuat di | Data seimbang | Data tidak seimbang |
| Sumbu Y | Sensitivitas (Recall) | Precision |
| Sumbu X | 1 - Spesifisitas | Recall |
Kurva Precision-Recall dapat divisualisasikan menggunakan berbagai
package di R, seperti PRROC atau precrec,
untuk menganalisis performa model dalam konteks ketidakseimbangan
data.
library(PRROC)
## Warning: package 'PRROC' was built under R version 4.3.3
## Loading required package: rlang
## Warning: package 'rlang' was built under R version 4.3.3
set.seed(57)
x1 <- rnorm(200)
x2 <- rbinom(200, 1, 0.5)
x3 <- rnorm(200)
lin_pred <- -1 + 1.3 * x1 - 0.5 * x2 + 0.7 * x3
p <- 1 / (1 + exp(-lin_pred))
y <- rbinom(200, 1, p)
data <- data.frame(y = y, x1, x2, x3)
model <- glm(y ~ x1 + x2 + x3, data = data, family = binomial)
prob <- predict(model, type = "response")
pr <- pr.curve(scores.class0 = prob[data$y == 1],
scores.class1 = prob[data$y == 0],
curve = TRUE)
plot(pr)
Pr Curve sangat informatif untuk aplikasi seperti deteksi penipuan atau diagnosis penyakit langka. PR Curve cocok digunakan saat: - Kelas positiv jauh lebih jarang dibanding kelas negatif - Tujuan aplikasi lebih mementingkan presisi terhadap kelas minoritas
Tujuan Pada sub-bab ini akan dijelaskan dan dihitung pseudo R-squared dalam regresi logistik: - \(R_{CoxandSnell}^2\) - \(R_{McFadden}^2\)
Simulasi Data
set.seed(57)
n <- 300
x1 <- rnorm(n)
x2 <- rbinom(n, 1, 0.5)
x3 <- rnorm(n)
lin_pred <- -1 + 1.4 * x1 - 0.4 * x2 + 0.8 * x3
p <- 1 / (1 + exp(-lin_pred))
y <- rbinom(n, 1, p)
data <- data.frame(y = as.factor(y), x1, x2, x3)
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_{CoxandSnell}^2=1-(\frac{L_0}{L_M})^{2/n} \] \[ R_{McFadden}^2=1-\frac{logL_M}{logL_0} \] Dengan: - \(L_0\): Likelihood model null (tanpa prediktor) - \(L_M\): Likelihood model penuh
Perhitungan Manual R-squared
logL0 <- logLik(model_null)
logLM <- logLik(model)
L0 <- exp(logL0)
LM <- exp(logLM)
n <- nobs(model)
cox_snell <- 1 - (L0 / LM)^(2 / n)
mcfadden <- 1 - (as.numeric(logLM) / as.numeric(logL0))
r2 <- data.frame(
R2_Cox_Snell = cox_snell,
R2_McFadden = mcfadden
)
r2
Perhitungan Otomatis dengan Package R Menggunakan pscl
library(pscl)
## Warning: package 'pscl' was built under R version 4.3.3
## Classes and Methods for R originally developed in the
## Political Science Computational Laboratory
## Department of Political Science
## Stanford University (2002-2015),
## by and under the direction of Simon Jackman.
## hurdle and zeroinfl functions by Achim Zeileis.
pR2(model)
## fitting null model for pseudo-r2
## llh llhNull G2 McFadden r2ML r2CU
## -124.6653690 -187.2993920 125.2680460 0.3344059 0.3413481 0.4786741
Menggunakan DescTools
if (!require(DescTools)) install.packages("DescTools"); library(DescTools)
PseudoR2(model, which = "all")
## McFadden McFaddenAdj CoxSnell Nagelkerke AldrichNelson
## 0.3344059 0.3130497 0.3413481 0.4786741 0.2945626
## VeallZimmermann Efron McKelveyZavoina Tjur AIC
## 0.5304650 0.3769812 0.5694989 0.3775846 257.3307379
## BIC logLik logLik0 G2
## 272.1458678 -124.6653690 -187.2993920 125.2680460
Interpretasi: - Nilai \(R^2\) mendekati 1 berarti model memiliki kekuatan prediktif yang baik - McFadden \(R^2\) > 0.2 sering dianggap sebagai model dengan kecocokan yang baik - Cox & Snell \(R^2\) lebih konservatif dan tidak pernah mencapai 1 penuh
Distribusi multinomial merupakan generalisasi dari distribusi binomial untuk kasus dengan lebih dari dua kategori. Jika suatu percobaan menghasilkan k kategori, dan kita ingin mengetahui probabilitas dari sejumlah kejadian dalam setiap kategori, maka distribusi multinomial bisa digunakan.
Jika variabel acak \(( X_1, X_2, ..., X_k)\) menyatakan jumlah kejadian pada masing-masing dari \(( k )\) kategori, maka distribusi peluangnya diberikan oleh:
\[ P(X_1 = x_1, ..., X_k = x_k) = \frac{n!}{x_1!x_2!...x_k!}p_1^{x_1}p_2^{x_2}...p_k^{x_k} \]
dengan syarat:
Sebuah survei dilakukan terhadap 10 responden yang diminta memilih satu dari tiga buah favorit: Apel (A), Jeruk (J), dan Pisang (P).
Hasil Survei: - Apel: 3 orang
- Jeruk: 5 orang
- Pisang: 2 orang
Probabilitas Teoritis Preferensi: - \(( p_A
= 0.3 )\)
- \(( p_J = 0.5 )\)
- \(( p_P = 0.2 )\)
Pertanyaan: Berapa peluang bahwa 10 orang akan terdiri dari 3 pemilih Apel, 5 pemilih Jeruk, dan 2 pemilih Pisang?
Perhitungan Manual di R
n <- 10
x <- c(3, 5, 2)
p <- c(0.3, 0.5, 0.2)
# 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.08505
Probabilitas bahwa 3 orang memilih Apel, 5 Jeruk, dan 2 Pisang (dengan proporsi preferensi 0.3, 0.5, dan 0.2) adalah 0.08505. Distribusi multinomial digunakan untuk menghitung peluang dalam percobaan dengan beberapa kategori hasil. Rumus dasarnya merupakan generalisasi dari binomial untuk lebih dari dua kategori.
Regresi logistik multinomial digunakan ketika kita ingin menganalisis hubungan antara variabel respon kategorik yang memiliki lebih dari dua kategori dengan satu atau lebih variabel prediktor.
Misalkan variabel respon \(( Y )\) memiliki \(( K )\) kategori, dan kita memilih salah satu kategori (biasanya kategori terakhir) sebagai acuan (baseline), maka model logit untuk kategori \(( j )\) terhadap kategori baseline \(( K )\) adalah:
\[ \log\left(\frac{P(Y = j)}{P(Y = K)}\right) = \beta_{0j} + \beta_{1j}x_1 + \cdots + \beta_{pj}x_p \]
untuk \(( j = 1, 2, ..., K-1 )\).
Model baseline-category logit membandingkan log odds (logit) dari setiap kategori terhadap satu kategori acuan (baseline). Bentuk umumnya adalah:
\[ \log\left(\frac{\pi_j}{\pi_c}\right), \quad j = 1, ..., c-1 \]
di mana:
Dengan demikian, akan terbentuk sebanyak \(( (c - 1) )\) fungsi logit.
Catatan: Dalam R, kategori referensi biasanya
otomatis ditetapkan sebagai kategori terakhir, namun bisa diubah
menggunakan fungsi relevel().
Model Regresi untuk Satu Prediktor
Jika hanya terdapat satu prediktor \(( x )\), maka bentuk persamaan regresi logit menjadi:
\[ \log\left(\frac{\pi_j}{\pi_c}\right) = \alpha_j + \beta_j x, \quad j = 1, ..., c-1 \]
Contoh Kasus: Respon dengan 3 Kategori
Misalnya variabel \(( Y )\) memiliki tiga kategori: \(( Y \in \{1, 2, 3\} )\), dan kategori ke-3 dijadikan acuan. Maka modelnya:
\[ \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 \]
Relasi Antar Kategori Untuk melihat hubungan antara kategori 1 dan 2, kita bisa menghitung logit relatifnya:
\[ \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 \]
Dengan ini, kita dapat membandingkan logit antara dua kategori non-baseline.
Implementasi di R Model ini dapat diimplementasikan
dengan fungsi multinom() dari package
nnet.
Parameter dalam model ini biasanya diestimasi menggunakan metode Maximum Likelihood, dengan pendekatan numerik seperti algoritma Newton-Raphson.
Fungsi log-likelihood untuk data dengan \(( n )\) observasi dan \(( k )\) kategori:
\[ \ell(\beta) = \sum_{i=1}^n \sum_{j=1}^k y_{ij} \log(\pi_{ij}) \]
dengan:
Sebuah restoran cepat saji ingin mengetahui faktor-faktor yang memengaruhi preferensi pelanggan terhadap jenis minuman yang dipesan, yaitu: Teh, Kopi, atau Jus.
Restoran melakukan survei terhadap 200 pelanggan dan mengumpulkan data berikut:
Tujuannya adalah: Mengetahui bagaimana usia, waktu kunjungan, dan frekuensi kunjungan memengaruhi pilihan jenis minuman pelanggan.
set.seed(57)
n <- 200
TimeofDay <- sample(c("Pagi", "Siang", "Malam"), n, replace = TRUE)
Age <- round(rnorm(n, mean = 30, sd = 8))
VisitFrequency <- round(pmax(rnorm(n, mean = 10, sd = 3), 0))
# Simulasikan Time of Day berdasarkan probabilitas berbeda per Jenis Minuman
Drink <- sapply(TimeofDay, function(time) {
if (time == "Pagi") {
sample(c("Teh", "Kopi", "Jus"), size = 1, prob = c(0.3, 0.4, 0.3))
} else if (time == "Siang") {
sample(c("Teh", "Kopi", "Jus"), size = 1, prob = c(0.2, 0.4, 0.4))
} else {
sample(c("Teh", "Kopi", "Jus"), size = 1, prob = c(0.4, 0.4, 0.2))
}
})
df <- data.frame(Drink = factor(Drink), Age, TimeofDay = factor(TimeofDay), VisitFrequency)
df$Drink <- relevel(df$Drink, ref = "Teh") # baseline
head(df)
library(nnet)
## Warning: package 'nnet' was built under R version 4.3.3
model_mnlogit <- multinom(Drink ~ Age + TimeofDay + VisitFrequency, data = df)
## # weights: 18 (10 variable)
## initial value 219.722458
## iter 10 value 209.542916
## final value 209.300869
## converged
summary(model_mnlogit)
## Call:
## multinom(formula = Drink ~ Age + TimeofDay + VisitFrequency,
## data = df)
##
## Coefficients:
## (Intercept) Age TimeofDayPagi TimeofDaySiang VisitFrequency
## Jus -1.26942303 0.0074754230 1.0855274 1.4624395 0.027840437
## Kopi -0.01732173 0.0004726801 0.7983403 0.8789951 -0.003789478
##
## Std. Errors:
## (Intercept) Age TimeofDayPagi TimeofDaySiang VisitFrequency
## Jus 0.9710145 0.02200367 0.4878910 0.4908274 0.06270705
## Kopi 0.8347962 0.01969198 0.4224498 0.4376284 0.05616052
##
## Residual Deviance: 418.6017
## AIC: 438.6017
z <- summary(model_mnlogit)$coefficients / summary(model_mnlogit)$standard.errors
pval <- 2 * (1 - pnorm(abs(z)))
round(pval, 4)
## (Intercept) Age TimeofDayPagi TimeofDaySiang VisitFrequency
## Jus 0.1911 0.7341 0.0261 0.0029 0.6571
## Kopi 0.9834 0.9808 0.0588 0.0446 0.9462
Interpretasi: Koefisien untuk kategori “Kopi” dan “Jus” dibandingkan dengan baseline “Teh”. Jika nilai p-value < 0.05 menunjukkan bahwa variabel tersebut signifikan memengaruhi preferensi minuman yang dipesan di restauran tersebut.
df$Predicted <- predict(model_mnlogit)
table(Predicted = df$Predicted, Actual = df$Drink)
## Actual
## Predicted Teh Jus Kopi
## Teh 26 12 25
## Jus 1 2 4
## Kopi 27 45 58
Model regresi logistik multinomial terbukti efektif dalam: - Mengkaji keterkaitan antara karakteristik pelanggan dan preferensi pemilihan jenis minuman
Mengidentifikasi faktor-faktor yang secara signifikan memengaruhi keputusan pemilihan
Memberikan kemampuan untuk memprediksi jenis minuman yang akan dipesan oleh pelanggan berdasarkan atribut pribadi mereka
Regresi logistik ordinal digunakan saat variabel respons memiliki sifat urutan, contohnya seperti tingkat kepuasan pelanggan: Tidak Puas, Cukup Puas, dan Sangat Puas.
Model ini berbeda dengan:
Model yang digunakan adalah Cumulative Logit Model dengan asumsi proportional odds. Model ini memodelkan probabilitas kumulatif hingga kategori tertentu:
\[ \log\left(\frac{P(Y \leq j)}{P(Y > j)}\right) = \alpha_j + \beta x \]
Penjelasan: - \(\alpha_j\): intercept khusus untuk batas kategori ke-\(j\) - \(\beta\): koefisien regresi (diasumsikan sama untuk semua batas kategori)
Jika ada \(k\) kategori, maka dibentuk \((k - 1)\) model logit kumulatif.
Koefisien \(\beta\) menunjukkan arah dan besar pengaruh dari variabel prediktor \(x\) terhadap peluang berada pada kategori lebih rendah atau sama.
Odds ratio (rasio peluang) dihitung sebagai:
\[ \text{OR} = e^{\beta} \] ## Contoh Kasus Sebuah kampus swasta melakukan survei terhadap mahasiswanya untuk mengevaluasi tingkat kepuasan terhadap fasilitas kampus. Tujuan dari survei ini adalah untuk mengetahui apakah jumlah fasilitas yang tersedia di kampus (seperti ruang belajar, Wi-Fi, kantin, laboratorium, area olahraga, dll.) berkaitan dengan kepuasan mahasiswa terhadap layanan dan kenyamanan lingkungan kampus.
Dalam penelitian ini, kampus dikelompokkan berdasarkan jumlah fasilitas yang tersedia menjadi tiga kategori: - Sedikit (1-3 fasilitas) - Sedang (4-6 fasilitas) - Banyak (7 fasilitas atau lebih)
Sementara itu, tingkat kepuasan mahasiswa dikategorikan menjadi tiga tingkat: - Tidak puas - Cukup - Puas
set.seed(57)
n <- 200
quality <- round(runif(n, 1, 10))
satisfaction <- cut(4 + 0.4*quality + rnorm(n),
breaks = c(-Inf, 5.5, 7.5, Inf),
labels = c("Tidak Puas", "Cukup", "Puas"),
ordered_result = TRUE)
df <- data.frame(satisfaction, quality)
head(df)
model_ord <- polr(satisfaction ~ quality, data = df, Hess = TRUE)
summary(model_ord)
## Call:
## polr(formula = satisfaction ~ quality, data = df, Hess = TRUE)
##
## Coefficients:
## Value Std. Error t value
## quality 0.5442 0.06487 8.389
##
## Intercepts:
## Value Std. Error t value
## Tidak Puas|Cukup 1.9201 0.3488 5.5053
## Cukup|Puas 4.6714 0.4916 9.5033
##
## Residual Deviance: 336.6941
## AIC: 342.6941
(ctable <- coef(summary(model_ord)))
## Value Std. Error t value
## quality 0.5442158 0.06486982 8.389352
## Tidak Puas|Cukup 1.9200523 0.34876696 5.505259
## Cukup|Puas 4.6713599 0.49155106 9.503305
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
## quality 0.5442158 0.06486982 8.389352 0
## Tidak Puas|Cukup 1.9200523 0.34876696 5.505259 0
## Cukup|Puas 4.6713599 0.49155106 9.503305 0
newdata <- data.frame(quality = 5:9)
predict(model_ord, newdata = newdata, type = "probs")
## Tidak Puas Cukup Puas
## 1 0.30980594 0.5656713 0.1245227
## 2 0.20664925 0.5964953 0.1968554
## 3 0.13130624 0.5717406 0.2969531
## 4 0.08064066 0.4981067 0.4212526
## 5 0.04843479 0.3951607 0.5564045
Model cumulative logit mengasumsikan efek prediktor sama untuk setiap cutof. Jika tidak, pertimbangkan model non-proportional odds seperti generalized ordinal model.
Selain model cumulative logit, terdapat beberapa pendekatan lain untuk menangani data ordinal:
Model alternatif ini bisa digunakan apabila asumsi proportional odds tidak terpenuhi pada cumulative logit.
polr() dari paket MASS.Regresi logistik ordinal yang paling umum digunakan adalah Cumulative Logit Model dengan asumsi proportional odds.
Asumsi ini juga dikenal sebagai asumsi paralelisme (parallel lines assumption), yang menyatakan bahwa: Koefisien regresi \(( \beta )\) tetap sama untuk setiap batas kumulatif dari kategori respons.
Rumus umum modelnya: \[ \log\left(\frac{P(Y \leq j)}{P(Y > j)}\right) = \alpha_j + \beta x,\quad j = 1, \ldots, c-1 \]
Ciri-cirinya: - Intercept \(( \alpha_j )\) berbeda untuk setiap batas kategori - Koefisien \(( \beta )\) tetap, menunjukkan efek variabel prediktor sama untuk semua kategori
Visualisasi asumsi paralelisme: - Jika asumsi ini terpenuhi, garis logit kumulatif terhadap \(( x )\) akan sejajar (slope sama, hanya beda posisi)
Jika Asumsi Tidak Terpenuhi: - Koefisien regresi bisa bervariasi di setiap batas kategori - Gunakan model yang lebih fleksibel, seperti: - Partial Proportional Odds Model - Generalized Ordinal Logistic Regression
Pengujian Asumsi Paralelisme
Untuk menguji apakah asumsi paralelisme terpenuhi, dapat digunakan:
brant di R)Contoh dengan package brant di R:
library(MASS)
library(brant)
## Warning: package 'brant' was built under R version 4.3.3
brant(model_ord)
## --------------------------------------------
## Test for X2 df probability
## --------------------------------------------
## Omnibus 0.46 1 0.5
## quality 0.46 1 0.5
## --------------------------------------------
##
## H0: Parallel Regression Assumption holds
Jika p-value < 0.05 maka asumsi nya tidak terpenuhi
Kesimpulan - Asumsi paralelisme penting untuk validitas model cumulative logit. - Menyederhanakan interpretasi karena efek prediktor konstan. - Jika tidak terpenuhi, gunakan model ordinal alternatif.
Analisis data kategorik berperan penting dalam statistik terapan karena banyak fenomena di dunia nyata menghasilkan data dalam bentuk kategori, seperti jenis kelamin, status pekerjaan, dan tingkat pendidikan.
Pendekatan statistik yang umum digunakan dalam analisis data kategorik meliputi:
Tabel Kontingensi
Digunakan untuk eksplorasi hubungan antar dua atau lebih variabel
kategorik dengan cara menyajikan frekuensi gabungan. Tabel ini bisa
dilengkapi dengan uji asosiasi seperti chi-square untuk mengetahui
hubungan antar variabel.
Model Log-linear
Digunakan untuk menganalisis struktur asosiasi dari tabel kontingensi
tanpa variabel dependen. Model ini cocok untuk memahami interaksi antar
variabel kategorik dan dapat diuji menggunakan likelihood ratio
test.
Model Regresi Logistik
Cocok digunakan saat ada variabel kategorik sebagai variabel dependen
(contoh: sakit/tidak). Model ini memodelkan peluang terjadinya suatu
kejadian (log odds) berdasarkan variabel independen dan dapat diperluas
menjadi regresi logistik multinomial atau ordinal.
Perbandingan
Kesimpulan
Pemilihan pendekatan tergantung pada tujuan analisis: apakah untuk deskripsi, eksplorasi struktur, atau prediksi. Kombinasi pendekatan sering digunakan untuk pemahaman yang lebih komprehensif terhadap data kategorik.
Tabel Kontingensi
Tabel kontingensi menyajikan jumlah frekuensi dari kombinasi kategori antar variabel.
Contoh:
table_data <- matrix(c(40, 15, 35, 60), nrow=2,
dimnames = list(Diet = c("Rendah Karbo", "Rendah Lemak"),
Hasil = c("Ya", "Tidak")))
table_data
## Hasil
## Diet Ya Tidak
## Rendah Karbo 40 35
## Rendah Lemak 15 60
Tabel kontingensi bersifat deskriptif dan tidak melibatkan pemodelan probabilitas.
Model Loglinear
Model loglinear memodelkan logaritma dari ekspektasi frekuensi sel dalam tabel kontingensi.
\[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \]
library(MASS)
loglm(~ Diet * Hasil, data = table_data)
## Call:
## loglm(formula = ~Diet * Hasil, data = table_data)
##
## Statistics:
## X^2 df P(> X^2)
## Likelihood Ratio 0 0 1
## Pearson 0 0 1
Model Regresi Logistik
Model regresi logistik biner:
\[ \log\left( \frac{p}{1 - p} \right) = \beta_0 + \beta_1 x \]
Digunakan jika ada variabel dependen kategorik (biasanya biner).
Bertujuan untuk memprediksi probabilitas suatu outcome.
Umumnya digunakan dalam studi observasional atau eksperimental. Contoh:
data_glm <- data.frame(
Hasil = c(1, 0, 1, 0),
Diet = factor(c("Rendah Karbo", "Rendah Karbo", "Rendah Lemak", "Rendah Lemak")),
Frek = c(40, 15, 35, 60)
)
model_logit <- glm(Hasil ~ Diet, weights = Frek, family = binomial, data = data_glm)
summary(model_logit)
##
## Call:
## glm(formula = Hasil ~ Diet, family = binomial, data = data_glm,
## weights = Frek)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.9808 0.3028 3.240 0.0012 **
## DietRendah Lemak -1.5198 0.3700 -4.108 4e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 207.94 on 3 degrees of freedom
## Residual deviance: 189.50 on 2 degrees of freedom
## AIC: 193.5
##
## Number of Fisher Scoring iterations: 4
| Aspek | Tabel Kontingensi | Model Loglinear | Regresi Logistik |
|---|---|---|---|
| Tujuan | Deskripsi frekuensi | Deteksi asosiasi | Prediksi probabilitas |
| Variabel Dependen | Tidak ada | Tidak ada (simetris) | 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
matrix(c(40, 15, 35, 60), nrow=2,
dimnames = list(Diet = c("Rendah Karbo", "Rendah Lemak"),
Hasil = c("Ya", "Tidak")))
## Hasil
## Diet Ya Tidak
## Rendah Karbo 40 35
## Rendah Lemak 15 60
Model log-linier untuk tabel I x J dapat dituliskan: \[ \log(\mu_{ij}) = \lambda + \lambda_i^{T} + \lambda_j^{R} + \lambda_{ij}^{TR} \]
Model saturated atau model penuh menyertakan seluruh efek utama dan interaksi:
# Data
library(MASS)
data <- matrix(c(45, 55, 65, 35), nrow=2, byrow=TRUE)
dimnames(data) <- list(Diet = c("Rendah Karbo", "Rendah Lemak"), Hasil = c("Ya", "Tidak"))
ftable(data)
## Hasil Ya Tidak
## Diet
## Rendah Karbo 45 55
## Rendah Lemak 65 35
Model saturated dapat dipasang dengan loglm dari package {MASS}:
model_saturated <- loglm(~ Diet * Hasil, data = data)
summary(model_saturated)
## Formula:
## ~Diet * Hasil
## attr(,"variables")
## list(Diet, Hasil)
## attr(,"factors")
## Diet Hasil Diet:Hasil
## Diet 1 0 1
## Hasil 0 1 1
## attr(,"term.labels")
## [1] "Diet" "Hasil" "Diet:Hasil"
## 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}) = \lambda + \lambda_i^{T} + \lambda_j^{R} \]
Model ini menguji hipotesis bahwa variabel X dan Y saling independen.
model_indep <- loglm(~ Diet + Hasil, data = data)
summary(model_indep)
## Formula:
## ~Diet + Hasil
## attr(,"variables")
## list(Diet, Hasil)
## attr(,"factors")
## Diet Hasil
## Diet 1 0
## Hasil 0 1
## attr(,"term.labels")
## [1] "Diet" "Hasil"
## 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 8.138435 1 0.004333667
## Pearson 8.080808 1 0.004473649
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.8197099
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:
## ~Diet + Hasil
## Model 2:
## ~Diet * Hasil
##
## Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1 8.138435 1
## Model 2 0.000000 0 8.138435 1 0.00433
## Saturated 0.000000 0 0.000000 0 1.00000
Model log-linear merupakan metode yang digunakan untuk menganalisis hubungan antara dua atau lebih variabel kategorik yang disusun dalam bentuk tabel kontingensi. Model ini bekerja dengan mengasumsikan bahwa logaritma dari nilai harapan frekuensi setiap sel (μᵢⱼ) dapat dinyatakan sebagai penjumlahan dari efek masing-masing variabel serta (jika diperlukan) efek interaksi antar variabel. Untuk tabel 2x2, bentuk modelnya adalah: \[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \] ## Perbedaan Utama Antara Model Log Linear dan Model Regresi Logistik - Model log-linear digunakan untuk memodelkan frekuensi (jumlah kejadian) dalam tabel kontingensi serta menguji hubungan antar variabel kategorik, tanpa menganggap adanya variabel respon dan prediktor.
Sistem Persamaan Model Log-Linear \[\begin{align*} \log(\mu_{11}) &= \lambda + \lambda^A_1 + \lambda^B_1 + \lambda^{AB}_{11} \\ \\ \log(\mu_{12}) &= \lambda + \lambda^A_1 + \lambda^B_2 + \lambda^{AB}_{12} \\ \\ \log(\mu_{21}) &= \lambda + \lambda^A_2 + \lambda^B_1 + \lambda^{AB}_{21} \\ \\ \log(\mu_{22}) &= \lambda + \lambda^A_2 + \lambda^B_2 + \lambda^{AB}_{22} \\ \end{align*}\] Constraint Sum-to-Zero \[\begin{align*} \lambda_1^A + \lambda_2^A = 0 \\ \\ \lambda_1^B + \lambda_2^B = 0 \\ \\ \lambda_{11}^{AB} + \lambda_{12}^{AB} + \lambda_{21}^{AB} + \lambda_{22}^{AB} = 0 \end{align*}\] Rumus Estimasi Parameter dengan Sum-to-Zero Constraint
\[ \lambda_1^A = \frac{1}{2} \left[ (\log \mu_{11} + \log \mu_{12}) - (\log \mu_{21} + \log \mu_{22}) \right] \] \[ \lambda_1^B = \frac{1}{2} \left[ (\log \mu_{11} + \log \mu_{21}) - (\log \mu_{12} + \log \mu_{22}) \right] \] \[ \lambda_{12}^{AB} = \frac{1}{4} \left[ \log \mu_{12} - \log \mu_{11} - \log \mu_{22} + \log \mu_{21} \right] \] ## Analisis Data Tabel Kontingensi 2x2
Seorang ilmuwan gizi melakukan penelitian untuk mengetahui apakah diet berpengaruh signifikan terhadap keberhasilan seseorang dalam menurunkan berat badan. Penelitian ini dilakukan terhadap 90 partisipan yang memiliki tujuan menurunkan berat badan dalam waktu 3 bulan.
Partisipan dibagi ke dalam dua kelompok berdasarkan kebiasaan mereka menjalankan program diet: - Melakukan Diet - Tidak Melakukan Diet
Setelah 3 bulan, dicatat apakah masing-masing partisipan berhasil menurunkan berat badan secara signifikan (misalnya, minimal 5% dari berat badan awal) atau tidak berhasil.
Data kemudian dikumpulkan dan disajikan dalam tabel kontingensi 2x2 berikut:
| Diet | Berhasil | Tidak |
|---|---|---|
| Ya | 40 | 15 |
| Tidak | 10 | 25 |
Model log-linear pada tabel 2x2: \[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \] dengan constraint sum-to-zero: \[ \sum_i \lambda_i^A = 0,\quad \sum_j \lambda_j^B = 0,\quad \sum_{i,j} \lambda_{ij}^{AB} = 0 \] ## Estimasi Parameter Model (Manual, Sum-to-zero) Misalkan: - A1 = Diet (Ya), A2 = Tidak - B1 = Berhasil, B2 = Tidak Observasi: - \(n_{11}\) = 40, \(n_{12}\) = 15 - \(n_{21}\) = 10, \(n_{22}\) = 25
\[\begin{align*} \log(\mu_{11}) &= \lambda + \lambda^A_1 + \lambda^B_1 + \lambda^{AB}_{11} \\ \\ \log(\mu_{12}) &= \lambda + \lambda^A_1 + \lambda^B_2 + \lambda^{AB}_{12} \\ \\ \log(\mu_{21}) &= \lambda + \lambda^A_2 + \lambda^B_1 + \lambda^{AB}_{21} \\ \\ \log(\mu_{22}) &= \lambda + \lambda^A_2 + \lambda^B_2 + \lambda^{AB}_{22} \\ \end{align*}\]
Constraint sum-to-zero:
\[\begin{align*} \lambda_1^A + \lambda_2^A = 0 \\ \\ \lambda_1^B + \lambda_2^B = 0 \\ \\ \lambda_{11}^{AB} + \lambda_{12}^{AB} + \lambda_{21}^{AB} + \lambda_{22}^{AB} = 0 \end{align*}\]
Langkah-langkah: 1. Hitung rata-rata log frekuensi sel: \[\begin{align*} \lambda &= \frac{1}{4} \sum_{i=1}^{2} \sum_{j=1}^{2} \log(n_{ij}) \\ &= \frac{1}{4} [\log(40) + \log(15) + \log(10) + \log(25)] \\ &= 2.979 \end{align*}\]
Efek Utama A (Merokok): \[\begin{align*} \lambda_1^A &= \frac{1}{2} [(\log(40) + \log(15)) - (\log(10) + \log(25))] \\ \\ &= \frac{1}{2} [(3.6889 + 2.7081) - (2.3026 + 3.2189)] \\ \\ &= \frac{1}{2} (6.3970 - 5.5215) \\ \\ &= \frac{1}{2} (0.8755) \\ \\ &= 0.4378 \\ \\ \lambda_2^A &= -0.4378 \\ \end{align*}\]
Efek Utama B (Status): \[\begin{align*} \\ \lambda_1^B &= \frac{1}{2} [(\log(40) + \log(10)) - (\log(15) + \log(25))] \\ \\ &= \frac{1}{2} [(3.6889 + 2.3026) - (2.7081 + 3.2189)] \\ \\ &= \frac{1}{2} (5.9915 - 5.9270) \\ \\ &= \frac{1}{2} (0.0645) \\ \\ &= 0.0323 \\ \\ \lambda_2^B &= -0.0323 \\ \end{align*}\]
Efek Interaksi: \[\begin{align*} \lambda_{11}^{AB} &= \frac{1}{4} [(\log(40) - \log(15)) - \log(10) + \log(25)] \\ \\ &= \frac{1}{4} [(3.6889 - 2.7081 - 2.3026 + 3.2189)] \\ \\ &= \frac{1}{4} (-1.3218 + 3.2189) \\ \\ &= \frac{1}{4} (1.8971) \\ \\ &= 0.4743 \\ \\ \lambda_{12}^{AB} &= -0.4743 \\ \\ \lambda_{21}^{AB} &= -0.4743 \\ \\ \lambda_{22}^{AB} &= 0.4743 \\ \end{align*}\]
Ringkasan Parameter
\[ OR = \frac{n_{11} n_{22}}{n_{12} n_{21}} = \frac{40 \times 25}{15 \times 10} = \frac{1000}{150} = 6.67 \] Log odds ratio: \[ \log(OR) = \log(6.67) = 1.8971 \] Standard error (SE):
\[ SE = \sqrt{ \frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}} } \\ \\ = \sqrt{ \frac{1}{40} + \frac{1}{15} + \frac{1}{10} + \frac{1}{25} }\\ \\ = \sqrt{0.025 + 0.067 + 0.1 + 0.04}\\ \\ = \sqrt{0.232}\\ \\ = 0.4817 \\ \]
95% Confidence Interval for log(OR):
\[\begin{align*} \\ \log(OR) \pm 1.96 \times SE &= 1.8971 \pm 1.96 \times 0.4817 \\ \\ &= 1.8971 \pm 0.9441 \\ \\ &= (1.8971 - 0.9441,\ 1.8971 + 0.9441) \\ \\ &= (0.9530,\ 2.8412) \\ \end{align*}\]
Back-transform to get CI for OR:
\[\begin{align*} \\ Lower &= \exp(0.9535) = 2.5935\\ \\ Upper &= \exp(2.8412) = 17.136 \\ \end{align*}\]
Jadi, OR=6.67 dengan taraf signifikansi 95% memiliki interval kepercayaan dalam rentang 2.5935 - 17.136
# Data 2x2
tabel <- matrix(c(40, 15, 10, 25), nrow = 2, byrow = TRUE)
colnames(tabel) <- c("Berhasil", "Tidak")
rownames(tabel) <- c("Ya", "Tidak")
tabel
## Berhasil Tidak
## Ya 40 15
## Tidak 10 25
#Ubah menjadi bentuk tabel untuk glm
data <- as.data.frame(as.table(tabel))
colnames(data) <- c("Diet", "Hasil", "Freq")
data
# Model tanpa interaksi
fit_no_inter <- glm(Freq ~ Diet + Hasil, family = poisson, data = data)
summary(fit_no_inter)
##
## Call:
## glm(formula = Freq ~ Diet + Hasil, family = poisson, data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.4195 0.1645 20.784 <2e-16 ***
## DietTidak -0.4520 0.2162 -2.090 0.0366 *
## HasilTidak -0.2231 0.2121 -1.052 0.2928
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 22.915 on 3 degrees of freedom
## Residual deviance: 17.319 on 1 degrees of freedom
## AIC: 42.628
##
## Number of Fisher Scoring iterations: 4
# Model dengan interaksi
fit_inter <- glm(Freq ~ Diet * Hasil, family = poisson, data = data)
summary(fit_inter)
##
## Call:
## glm(formula = Freq ~ Diet * Hasil, family = poisson, data = data)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.6889 0.1581 23.331 < 2e-16 ***
## DietTidak -1.3863 0.3536 -3.921 8.82e-05 ***
## HasilTidak -0.9808 0.3028 -3.240 0.0012 **
## DietTidak:HasilTidak 1.8971 0.4813 3.942 8.10e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 2.2915e+01 on 3 degrees of freedom
## Residual deviance: 4.4409e-16 on 0 degrees of freedom
## AIC: 27.309
##
## Number of Fisher Scoring iterations: 3
Nilai \(\log(6)=1.8971\) sama dengan efek interaksi pada output R
Suatu survei dilakukan untuk mengetahui hubungan antara Jenis Kelamin (Laki-laki/Perempuan) dan Kategori BMI (Kurus/Normal/Gemuk):
| Rendah | Sedang | Tinggi | |
|---|---|---|---|
| Pekerja | 10 | 25 | 15 |
| Pengangguran | 20 | 18 | 5 |
Bentuk umum model log-linear untuk tabel 2x3 (dengan sum-to-zero constraint): \[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \] dengan: - \(\mu_{ij}\): ekspektasi frekuensi pada baris ke-i, kolom ke-j - A: Status Pekerjaan (i=1:Pekerja, i=2:Pengangguran) - B: Tingkat Kepuasan Hidup (j=1:Rendah, j=2:Sedang, j=3:Tinggi) - Constraint: \[ \sum_i \lambda_i^A = 0,\quad \sum_j \lambda_j^B = 0,\quad \sum_{i} \lambda_{ij}^{AB} = 0, dan \sum_{j}\lambda_{ij}^{AB}=0 \] Secara eksplisit:
\[\begin{align*} \\ \log(\mu_{ij}) &= \lambda \\ \\ &+ \lambda^A_1(Pekerja), \lambda_2^A(Pengangguran)\\ &+ \lambda^B_1(Rendah), \lambda^B_2(Sedang), \lambda^B_3(Tiggi)\\ \\ &+ \lambda^{AB}_{ij}(interaksi jika ada) \\ \end{align*}\]
Langkah-langkah perhitungan:
Hitung rata-rata log frekuensi sel: \[\begin{align*} \\ \lambda &= \frac{1}{6} \sum_{i=1}^{2} \sum_{j=1}^{3} \log(n_{ij}) \\ \\ &= \frac{1}{6} [\log(10) + \log(25) + \log(15) + \log(20) + \log(18) + \log(5)] \\ \\ &= \frac{1}{6} [2.3026 + 3.2189 + 2.7081 + 2.9957 + 2.8904 + 0.6931] \\ \\ &= \frac{1}{6} (14.8088) \\ \\ &= 2.4681 \\ \end{align*}\]
Efek utama A (Jenis Kelamin): \[\begin{align*} \\ \lambda_1^A &= \frac{1}{3} \left[ (\log 10 + \log 25 + \log 15) - (\log 20 + \log 18 + \log 5) \right] \\ \\ &= \frac{1}{3} [ (2.3026 + 3.2189 + 2.7081) - (2.9957 + 2.8904 + 0.6931) ] \\ \\ &= \frac{1}{3} (8.2296 - 6.5792) \\ \\ &= \frac{1}{3} (1.6504) \\ \\ &= 0.5501 \\ \\ \lambda_2^A &= -\lambda_1^A = -0.5501 \\ \end{align*}\]
Efek utama B (Kategori BMI): \[\begin{align*} \\ \lambda_1^B &= \frac{1}{2} \left[ (\log 10 + \log 20) - (\log 25 + \log 18 + \log 15 + \log 5)/2 \right] \\ \\ &= \frac{1}{2} (2.3026 + 2.9957 - 5.2134) \\ \\ &= \frac{1}{2} (0.0849) \\ \\ &= 0.0425 \\ \\ \lambda_2^B &= \frac{1}{2} (3.2189 + 2.8904 - 4.8079) = 0.6498 \\ \\ \lambda_3^B &= -\lambda_1^B - \lambda_2^B = -0.7347 \\ \end{align*}\]
Efek Interaksi: \[\begin{align*} \\ \lambda_{11}^{AB} &= \log(10) - \lambda - \lambda_1^A - \lambda_1^B \\ \\ &= 2.3026 - 2.4681 - 0.551 - 0.0425 = -0.759 \\ \\ \lambda_{12}^{AB} &= \log(25) - \lambda - \lambda_1^A - \lambda_2^B \\ \\ &= 3.2189 - 2.4681 - 0.551 - 0.6498 = -0.45 \\ \\ \lambda_{13}^{AB} &= \log(15) - \lambda - \lambda_1^A - \lambda_3^B \\ \\ &= 2.7081 - 2.4681 - 0.551 + 0.7347 = 0.4237 \\ \\ \lambda_{21}^{AB} &= \log(20) - \lambda - \lambda_2^A - \lambda_1^B \\ \\ &= 2.9957 - 2.4681 + 0.551 - 0.0425 = 1.0361 \\ \\ \lambda_{22}^{AB} &= \log(18) - \lambda - \lambda_2^A - \lambda_2^B \\ \\ &= 2.8904 - 2.4681 + 0.551 - 0.6494 = 0.3239 \\ \\ \lambda_{23}^{AB} &= \log(5) - \lambda - \lambda_2^A - \lambda_3^B \\ \\ &= 0.6931 - 2.4681 + 0.551 + 0.7347 = -0.4893 \\ \end{align*}\]
Hitung Odds Ratio dan Interval Kepercayaan
\[ SE = \sqrt{\frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}}}\\ \\ = \sqrt{\frac{1}{10} + \frac{1}{18} + \frac{1}{20} + \frac{1}{25}}\\ \\ = \sqrt{0.1 + 0.0556 + 0.05 + 0.04}\\ \\ = \sqrt{0.2456}\\ \\ = 0.4956 \\ \]
Jadi, OR=0.36 dengan taraf signifikansi 95% memiliki interval kepercayaan dengan rentang 0.1363 - 0.9509
# Membuat data frame dari tabel
tabel2x3 <- matrix(c(10, 25, 15, 20, 18, 5), nrow = 2, byrow = TRUE)
colnames(tabel2x3) <- c("Rendah", "Sedang", "Tinggi")
rownames(tabel2x3) <- c("Pekerja", "Pengangguran")
tabel2x3
## Rendah Sedang Tinggi
## Pekerja 10 25 15
## Pengangguran 20 18 5
# Ubah menjadi data.frame untuk glm
data2x3 <- as.data.frame(as.table(tabel2x3))
colnames(data2x3) <- c("Status", "TingkatKepuasan", "Freq")
data2x3
# Model log-linear tanpa interaksi (asumsi independen)
fit_no_inter23 <- glm(Freq ~ Status + TingkatKepuasan, family = poisson, data = data2x3)
summary(fit_no_inter23)
##
## Call:
## glm(formula = Freq ~ Status + TingkatKepuasan, family = poisson,
## data = data2x3)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 2.7806 0.2064 13.475 <2e-16 ***
## StatusPengangguran -0.1508 0.2080 -0.725 0.468
## TingkatKepuasanSedang 0.3600 0.2379 1.513 0.130
## TingkatKepuasanTinggi -0.4055 0.2887 -1.405 0.160
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 18.4178 on 5 degrees of freedom
## Residual deviance: 9.2477 on 2 degrees of freedom
## AIC: 44.085
##
## Number of Fisher Scoring iterations: 4
# Model log-linear dengan interaksi (untuk cek asosiasi)
fit_inter23 <- glm(Freq ~ Status * TingkatKepuasan, family = poisson, data = data2x3)
summary(fit_inter23)
##
## Call:
## glm(formula = Freq ~ Status * TingkatKepuasan, family = poisson,
## data = data2x3)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 2.3026 0.3162 7.281 3.3e-13
## StatusPengangguran 0.6931 0.3873 1.790 0.07350
## TingkatKepuasanSedang 0.9163 0.3742 2.449 0.01433
## TingkatKepuasanTinggi 0.4055 0.4082 0.993 0.32062
## StatusPengangguran:TingkatKepuasanSedang -1.0217 0.4955 -2.062 0.03924
## StatusPengangguran:TingkatKepuasanTinggi -1.7918 0.6455 -2.776 0.00551
##
## (Intercept) ***
## StatusPengangguran .
## TingkatKepuasanSedang *
## TingkatKepuasanTinggi
## StatusPengangguran:TingkatKepuasanSedang *
## StatusPengangguran:TingkatKepuasanTinggi **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 1.8418e+01 on 5 degrees of freedom
## Residual deviance: -2.8866e-15 on 0 degrees of freedom
## AIC: 38.838
##
## Number of Fisher Scoring iterations: 3
Interpretasi: Model tanpa interaksi: - Jika deviance tidak signifikan, maka Status Pekerjaan dan Tingkat Kepuasan Hidup independen.
Model dengan interaksi:
Contoh interpretasi hasil (misal): - Jika koefisien StatusPekerjaanPengangguran negatif: proporsi Pengangguran pada kategori referensi lebih kecil dibanding Pekerja. - Jika koefisien TingkatKepuasan_Sedang positif: kemungkinan seseorang dengan tingkat kepuasan Sedang lebih tinggi daripada Rendah (pada Pekerja). - Jika model interaksi signifikan, pola distribusi Tingkat Kepuasan berbeda antara Pekerja dan Pengangguran.
Pada bagian sebelumnya, telah dijelaskan bahwa salah satu tujuan utama dari penggunaan model log-linear adalah untuk memperkirakan parameter-parameter yang menggambarkan hubungan antar variabel kategorik.
Dalam materi kali ini, kita akan mempelajari model log-linear yang lebih kompleks, yaitu yang digunakan untuk tabel kontingensi tiga arah. Model ini melibatkan tiga variabel kategorik, sehingga kemungkinan munculnya interaksi antar variabel menjadi lebih banyak. Dalam hal ini, interaksi paling kompleks yang dapat dianalisis adalah interaksi tiga arah, yaitu interaksi yang melibatkan ketiga variabel sekaligus.
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).
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.
Model Conditional
Conditional pada X: \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \] Memuat interaksi X dengan Y dan X dengan Z.
Conditional pada Y: \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \] Memuat interaksi Y dengan X dan Y dengan Z.
Conditional pada Z: \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \] Memuat interaksi Z dengan X dan Z dengan Y.
Model Joint Independence
Independensi antara X & Y: \(\log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij}\)
Independensi antara X & Z: \(\log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XZ}_{ik}\)
Independensi antara Y & Z: \(\log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{YZ}_{jk}\)
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:
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.
Deskripsi Data
Sebuah tim peneliti psikologi dari universitas ternama ingin memahami bagaimana faktor-faktor seperti jenis kelamin dan frekuensi olahraga berkaitan dengan tingkat stres pada individu. Stres merupakan salah satu isu kesehatan mental yang umum, terutama di kalangan mahasiswa dan pekerja muda, sehingga penting untuk memahami faktor-faktor yang mungkin berkontribusi terhadap tingginya atau rendahnya tingkat stres seseorang.
Peneliti mengklasifikasikan responden ke dalam tiga dimensi utama: Tingkat Stres: - Tinggi - Sedang - Rendah
Jenis Kelamin: - Laki-laki - Perempuan
Frekuensi Olahraga: - Jarang (kurang dari 1x per minggu) - Sedang (1-2x per minggu) - Rutin (lebih dari 3x per minggu)
Survei dilakukan kepada 462 responden dari berbagai latar belakang. Data dikumpulkan dan disusun dalam format tiga dimensi, dan berikut adalah salah satu cara menyajikannya, yaitu dalam bentuk tabel silang per jenis kelamin.
Tabel data:
| Tingkat Stres | Jenis Kelamin | Jarang | Rutin | Total |
|---|---|---|---|---|
| Tinggi | Laki-laki | 52 | 20 | 72 |
| Perempuan | 60 | 16 | 76 | |
| Total | 112 | 36 | 148 | |
| Sedang | Laki-laki | 30 | 38 | 68 |
| Perempuan | 36 | 36 | 72 | |
| Total | 66 | 74 | 140 | |
| Rendah | Laki-laki | 18 | 62 | 80 |
| Perempuan | 14 | 80 | 94 | |
| Total | 32 | 142 | 174 |
library("epitools")
library("DescTools")
library("lawstat")
## Warning: package 'lawstat' was built under R version 4.3.3
# Input data sesuai tabel praktikum
z.stres <- factor(rep(c("1tinggi", "2sedang", "3rendah"), each = 4))
x.sex <- factor(rep(c("1L", "2P"), each = 2, times = 3))
y.freq <- factor(rep(c("1jarang", "2rutin"), times = 6))
counts <- c(52, 20, 60, 16, 30, 38, 36, 36, 18, 62, 14, 80)
data3 <- data.frame(
TingkatStres = z.stres,
JenisKelamin = x.sex,
FrekuensiOlahraga = y.freq,
Frekuensi = counts
)
data3
Membentuk Tabel Kontingensi 3 Arah
table3d <- xtabs(Frekuensi ~ TingkatStres + JenisKelamin + FrekuensiOlahraga, data = data3)
ftable(table3d)
## FrekuensiOlahraga 1jarang 2rutin
## TingkatStres JenisKelamin
## 1tinggi 1L 52 20
## 2P 60 16
## 2sedang 1L 30 38
## 2P 36 36
## 3rendah 1L 18 62
## 2P 14 80
Analisis Log-Linear: Tahap Pemodelan Kita akan memodelkan tabel ini menggunakan beberapa model log-linear dan membandingkan kecocokan model (parsimonious model):
Penentuan Kategori Referensi
##=============================##
# Penentuan kategori reference
##=============================##
x.sex <- relevel(x.sex, ref = "2P")
y.freq <- relevel(y.freq, ref = "2rutin")
z.stres <- relevel(z.stres, ref = "3rendah")
Model Saturated memasukkan semua interaksi hingga tiga arah (tiga kategori): \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} + \lambda^{XYZ}_{ijk} \]
# Model saturated
model_saturated <- glm(counts ~ x.sex + y.freq + z.stres +
x.sex*y.freq + x.sex*z.stres + y.freq*z.stres +
x.sex*y.freq*z.stres,
family = poisson(link = "log"))
summary(model_saturated)
##
## Call:
## glm(formula = counts ~ x.sex + y.freq + z.stres + x.sex * y.freq +
## x.sex * z.stres + y.freq * z.stres + x.sex * y.freq * z.stres,
## family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.3820 0.1118 39.194 < 2e-16 ***
## x.sex1L -0.2549 0.1692 -1.506 0.132
## y.freq1jarang -1.7430 0.2897 -6.016 1.78e-09 ***
## z.stres1tinggi -1.6094 0.2739 -5.877 4.18e-09 ***
## z.stres2sedang -0.7985 0.2007 -3.979 6.93e-05 ***
## x.sex1L:y.freq1jarang 0.5062 0.3945 1.283 0.199
## x.sex1L:z.stres1tinggi 0.4780 0.3757 1.272 0.203
## x.sex1L:z.stres2sedang 0.3090 0.2876 1.074 0.283
## y.freq1jarang:z.stres1tinggi 3.0647 0.4039 7.589 3.23e-14 ***
## y.freq1jarang:z.stres2sedang 1.7430 0.3735 4.667 3.06e-06 ***
## x.sex1L:y.freq1jarang:z.stres1tinggi -0.8725 0.5514 -1.582 0.114
## x.sex1L:y.freq1jarang:z.stres2sedang -0.7426 0.5204 -1.427 0.154
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 1.2499e+02 on 11 degrees of freedom
## Residual deviance: 2.5757e-14 on 0 degrees of freedom
## AIC: 88.183
##
## Number of Fisher Scoring iterations: 3
exp(model_saturated$coefficients)
## (Intercept) x.sex1L
## 80.0000000 0.7750000
## y.freq1jarang z.stres1tinggi
## 0.1750000 0.2000000
## z.stres2sedang x.sex1L:y.freq1jarang
## 0.4500000 1.6589862
## x.sex1L:z.stres1tinggi x.sex1L:z.stres2sedang
## 1.6129032 1.3620072
## y.freq1jarang:z.stres1tinggi y.freq1jarang:z.stres2sedang
## 21.4285714 5.7142857
## x.sex1L:y.freq1jarang:z.stres1tinggi x.sex1L:y.freq1jarang:z.stres2sedang
## 0.4179259 0.4758772
Model yang digunakan adalah model log-linear saturated dengan semua efek utama, interaksi dua arah, dan interaksi tiga arah. Model ini memodelkan hubungan antara jenis kelamin (x.sex), frekuensi olahraga (y.freq), dan tingkat stress (z.stres) terhadap frekuensi responden. Dengan hasil estimasi koefisien nya sebagai berikut:
| Parameter | Estimate | Std. Error | z value | Pr(>) |
|---|---|---|---|---|
| (Intercept) | 4.3820 | 0.1118 | 39.194 | <2e-16 *** |
| x.sex1L | -0.2549 | 0.1692 | -1.506 | 0.132 |
| y.freq1jarang | -1.7430 | 0.2897 | -6.016 | 1.78e-09 *** |
| z.stres1tinggi | -1.6094 | 0.2739 | -5.877 | 4.18e-09 *** |
| z.stres2sedang | -0.7985 | 0.2007 | -3.979 | 6.93e-05 *** |
| x.sex1L:y.freq1jarang | 0.5062 | 0.3945 | 1.283 | 0.199 |
| x.sex1L:z.stres1tinggi | 0.4780 | 0.3757 | 1.272 | 0.203 |
| x.sex1L:z.stres2sedang | 0.3090 | 0.2876 | 1.074 | 0.283 |
| y.freq1jarang:z.stres1tinggi | 3.0647 | 0.4309 | 7.108 | 3.23e-14 *** |
| y.freq1jarang:z.stres2sedang | 1.7430 | 0.3735 | 4.667 | 3.06e-06 *** |
| x.sex1L:y.freq1jarang:z.stres1tinggi | -0.8725 | 0.5514 | -1.582 | 0.114 |
| x.sex1L:y.freq1jarang:z.stres2sedang | -0.7426 | 0.5204 | -1.427 | 0.154 |
Interpretasi Koefisien:
(Intercept):
Rata-rata log jumlah kasus untuk kategori referensi (Perempuan,
Frekuensi Sering, Stres Rendah) adalah 4.382. Ini menunjukkan baseline
expected count sebesar exp(4.382) ≈ 80.
x.sex1L (Laki-laki):
Laki-laki memiliki expected count sekitar
exp(-0.2549) ≈ 0.775 kali perempuan, tetapi perbedaannya
tidak signifikan (p = 0.132).
y.freq1jarang (Frekuensi olahraga jarang):
Responden dengan frekuensi olahraga yang jarang memiliki expected count
sekitar exp(-1.7430) ≈ 0.175 dari mereka yang sering
(signifikan, p < 0.001).
z.stres1tinggi (Stres tinggi):
Responden dengan tingkat stres tinggi memiliki expected count sekitar
exp(-1.6094) ≈ 0.20 dari mereka yang stres rendah
(signifikan, p < 0.001).
z.stres2sedang (Stres sedang):
Responden dengan stres sedang memiliki expected count sekitar
exp(-0.7985) ≈ 0.45 dari yang stres rendah (signifikan, p
< 0.001).
x.sex1L:y.freq1jarang:
Interaksi ini menunjukkan bahwa efek frekuensi olahraga jarang pada
laki-laki sedikit lebih tinggi (exp(0.5062) ≈ 1.66), tetapi
tidak signifikan (p = 0.199).
x.sex1L:z.stres1tinggi:
Interaksi ini menunjukkan bahwa stres tinggi pada laki-laki meningkatkan
expected count (exp(0.4780) ≈ 1.61), tapi tidak signifikan
(p = 0.203).
x.sex1L:z.stres2sedang:
Interaksi stres sedang pada laki-laki (exp(0.3090) ≈ 1.36)
juga tidak signifikan (p = 0.283).
y.freq1jarang:z.stres1tinggi:
Interaksi antara frekuensi olahraga jarang dan stres tinggi meningkatkan
expected count sebesar exp(3.0647) ≈ 21.4, sangat
signifikan (p < 0.001).
y.freq1jarang:z.stres2sedang:
Interaksi antara olahraga jarang dan stres sedang menghasilkan expected
count exp(1.7430) ≈ 5.7 kali lipat, signifikan (p <
0.001).
x.sex1L:y.freq1jarang:z.stres1tinggi:
Interaksi tiga arah ini menurunkan expected count
exp(-0.8725) ≈ 0.42, tetapi tidak signifikan (p =
0.114).
x.sex1L:y.freq1jarang:z.stres2sedang:
Interaksi tiga arah ini juga menurunkan expected count
exp(-0.7426) ≈ 0.48, tidak signifikan (p = 0.154).
Goodness-of-Fit
Residual deviance ≈ 0:
Menandakan bahwa model saturated sangat cocok dengan data (seluruh
variasi data dijelaskan model).
AIC = 88.18:
Dapat digunakan untuk perbandingan dengan model yang lebih sederhana.
Semakin kecil AIC, semakin baik model.
Model ini sangat fit terhadap data (residual deviance mendekati nol), tetapi tidak semua parameter signifikan.
Kesimpulan:
Model ini sangat fit terhadap data (residual deviance mendekati nol), tetapi tidak semua parameter signifikan.
Efek utama yang paling signifikan adalah:
Frekuensi olahraga jarang (expected count turun drastis dibandingkan yang sering)
Tingkat stres tinggi dan sedang (expected count turun signifikan dibandingkan stres rendah)
Interaksi antara olahraga jarang dan stres tinggi/sedang meningkatkan expected count secara signifikan
Tidak ditemukan bukti kuat untuk interaksi tiga arah yang signifikan
Model yang lebih sederhana (tanpa interaksi tiga arah) perlu dipertimbangkan untuk model final yang lebih parsimonious.
Catatan Interpretasi:
Nilai exp(coef) menyatakan rasio ekspektasi (expected count ratio) dibandingkan baseline.
Efek positif -> menaikkan expected count; Efek negatif -> menurunkan expected count.
Koefisien signifikan pada p-value < 0.05
Model log-linear homogenous memasukkan semua efek utama dan semua interaksi dua arah, tanpa interaksi tiga arah. Secara matematis, model ini dapat dituliskan sebagai berikut: \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
# Homogenous Model
model_homogenous <- glm(counts ~ x.sex + y.freq + z.stres +
x.sex*y.freq + x.sex*z.stres + y.freq*z.stres,
family = poisson(link = "log"))
summary(model_homogenous)
##
## Call:
## glm(formula = counts ~ x.sex + y.freq + z.stres + x.sex * y.freq +
## x.sex * z.stres + y.freq * z.stres, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.33523 0.11098 39.064 < 2e-16 ***
## x.sex1L -0.15079 0.15706 -0.960 0.337017
## y.freq1jarang -1.46409 0.21804 -6.715 1.88e-11 ***
## z.stres1tinggi -1.43945 0.22433 -6.417 1.39e-10 ***
## z.stres2sedang -0.70931 0.18283 -3.880 0.000105 ***
## x.sex1L:y.freq1jarang -0.05711 0.21327 -0.268 0.788866
## x.sex1L:z.stres1tinggi 0.13993 0.25524 0.548 0.583531
## x.sex1L:z.stres2sedang 0.12054 0.23561 0.512 0.608923
## y.freq1jarang:z.stres1tinggi 2.62706 0.27403 9.587 < 2e-16 ***
## y.freq1jarang:z.stres2sedang 1.37740 0.25888 5.321 1.03e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 124.9856 on 11 degrees of freedom
## Residual deviance: 2.9849 on 2 degrees of freedom
## AIC: 87.168
##
## Number of Fisher Scoring iterations: 4
Pengujian ini menggunakan residual deviance dari kedua model (saturated dan homogenous).
Langkah-Langkah Pengujian: 1. Hipotesis: - H0: Tidak ada interaksi tiga arah (model homogenous sudah cukup) - H1: Ada interaksi tiga arah (model saturated diperlukan)
# Deviance antar model
Deviance.model <- model_homogenous$deviance - model_saturated$deviance
Deviance.model
## [1] 2.984885
# Derajat bebas = db model homogenous - db model saturated
derajat.bebas <- (model_homogenous$df.residual - model_saturated$df.residual)
derajat.bebas
## [1] 2
chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel
## [1] 5.991465
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan
## [1] "Terima H0"
Interpretasi Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara jenis kelamin, frekuensi olahraga, dan tingkat stress.
Catatan: - Model pengurang adalah model yang lebih lengkap (lebih banyak parameter, df lebih kecil), yaitu model saturated - Derajat bebas dihitung dari selisih derajat bebas model homogenous dan saturated. - Keputusan berdasarkan perbandingan deviance model dengan chi-square tabel.
Rangkuman
Pengujian ada tidaknya interaksi tiga arah (Saturated Model vs Homogenous Model)
Hipotesis
H0: \(\lambda_{ijk}^{XYZ}=0\) (Tidak ada interaksi tiga arah; model yang terbentuk adalah model homogenous) -
H0: \(\lambda_{ijk}^{XYZ}\neq0\) (Tidak ada interaksi tiga arah; model yang terbentuk adalah model homogenous)
Tingkat Signifikansi
\(\alpha=5\%\)
Statistik Uji
\(\Delta\)Deviance = Deviance model homogenous - Deviance model saturated = 2.9849 - 0.00 = 2.9849
\(db\) = \(db\) model homogenous - \(db\) model saturated = 2 - 0 = 2
Daerah Penolakan
Tolak H0 jika \(\Delta\)Deviance \(\>>\) \(\chi_{0.05,db}^2 = \chi_{0.05,2}^2 = 5.991\)
Keputusan
Karena 2.9849 < 5.991, maka terima H0
Interpretasi
Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara jenis kelamin, frekuensi olahraga, dan tingkat stress.
Model log-linear conditional pada X memasukkan efek utama dan interaksi dua arah antara X dengan Y dan X dengan Z, tanpa interaksi antara Y dengan Z maupun interaksi tiga arah. \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} \]
# Conditional Association on X
model_conditional_X <- glm(counts ~ x.sex + y.freq + z.stres +
x.sex*y.freq + x.sex*z.stres,
family = poisson(link = "log"))
summary(model_conditional_X)
##
## Call:
## glm(formula = counts ~ x.sex + y.freq + z.stres + x.sex * y.freq +
## x.sex * z.stres, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.937e+00 1.187e-01 33.178 <2e-16 ***
## x.sex1L -1.613e-01 1.743e-01 -0.925 0.3548
## y.freq1jarang -1.823e-01 1.291e-01 -1.412 0.1579
## z.stres1tinggi -2.126e-01 1.543e-01 -1.378 0.1682
## z.stres2sedang -2.666e-01 1.566e-01 -1.702 0.0887 .
## x.sex1L:y.freq1jarang -3.517e-15 1.871e-01 0.000 1.0000
## x.sex1L:z.stres1tinggi 1.072e-01 2.240e-01 0.479 0.6323
## x.sex1L:z.stres2sedang 1.041e-01 2.274e-01 0.458 0.6471
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 124.99 on 11 degrees of freedom
## Residual deviance: 115.77 on 4 degrees of freedom
## AIC: 195.95
##
## Number of Fisher Scoring iterations: 5
Pengujian Ada Tidaknya Interaksi Antara Y dan Z (Homogenous Model vs Conditional Association on X)
Hipotesis
H0: \(\lambda_{jk}^{YZ}=0\) (Tidak ada interaksi antara frekuensi olahraga (Y) dan tingkat stress (Z)) -
H0: \(\lambda_{jk}^{YZ}\neq0\) (Ada interaksi antara frekuensi olahraga (Y) dan tingkat stress (Z))
Tingkat Signifikansi
\(\alpha=5\%\)
Statistik Uji
\(\Delta\)Deviance = Deviance model conditional on X - Deviance model homogenous = 115.77 - 2.9849 = 112.7855
\(db\) = \(db\) model homogenous - \(db\) model saturated = 4 - 2 = 2
Daerah Penolakan
Tolak H0 jika \(\Delta\)Deviance \(\>>\) \(\chi_{0.05,db}^2 = \chi_{0.05,2}^2 = 5.991\)
Keputusan
Karena 112.7855 < 5.991, maka tolak H0
Interpretasi
Dengan taraf nyata 5%, sudah cukup bukti untuk menolak H0 atau dapat dikatakan bahwa ada interaksi antara frekuensi olahraga dan tingkat stress. Dengan kata lain, model yang terbentuk adalah model dengan parameter \(\lambda_{jk}^{YZ}\)
# Pengujian hipotesis
# Deviance of Model
Deviance.model <- model_conditional_X$deviance - model_homogenous$deviance # model_conditional_X: conditional on X, model_homogenous: homogenous
Deviance.model
## [1] 112.7855
# Chi Square tabel dengan alpha = 0.05
derajat.bebas <- (4 - 2)
derajat.bebas
## [1] 2
chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 5.991465
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Tolak"
Interpretasi Karena nilai Deviance.model = 112.7855 lebih kecil dari nilai kritis chi-square tabel = 5.99 (dengan df = 2, alpha = 0.05), maka keputusan uji adalah “Tolak H0”.
Pada taraf nyata 5%, sudah cukup bukti untuk menolak H0, atau dengan kata lain ada interaksi antara frekuesi olahraga (Y) dan tingkat stress (Z). Model yang terbentuk adalah model homogenous, sehingga interaksi Y*Z signifikan secara statistik.
Model log-linear conditional pada Y memasukkan efek utama dan interaksi dua arah antara X dengan Y dan Y dengan Z, tanpa interaksi antara X dengan Z maupun interaksi tiga arah. \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{YZ}_{jk} \]
# Conditional Association on Y
model_conditional_Y <- glm(counts ~ x.sex + y.freq + z.stres +
x.sex*y.freq + y.freq*z.stres,
family = poisson(link = "log"))
summary(model_conditional_Y)
##
## Call:
## glm(formula = counts ~ x.sex + y.freq + z.stres + x.sex * y.freq +
## y.freq * z.stres, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.309e+00 1.032e-01 41.757 < 2e-16 ***
## x.sex1L -9.531e-02 1.261e-01 -0.756 0.45
## y.freq1jarang -1.490e+00 2.150e-01 -6.930 4.20e-12 ***
## z.stres1tinggi -1.372e+00 1.866e-01 -7.354 1.92e-13 ***
## z.stres2sedang -6.518e-01 1.434e-01 -4.546 5.47e-06 ***
## x.sex1L:y.freq1jarang 2.851e-15 1.871e-01 0.000 1.00
## y.freq1jarang:z.stres1tinggi 2.625e+00 2.739e-01 9.585 < 2e-16 ***
## y.freq1jarang:z.stres2sedang 1.376e+00 2.588e-01 5.316 1.06e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 124.9856 on 11 degrees of freedom
## Residual deviance: 3.3587 on 4 degrees of freedom
## AIC: 83.541
##
## Number of Fisher Scoring iterations: 4
Pengujian Ada Tidaknya Interaksi Antara X dan Z (Homogenous Model vs Conditional Association on Y)
Hipotesis
H0: \(\lambda_{ik}^{XZ}=0\) (Tidak ada interaksi antara jenis kelamin (X) dan tingkat stress (Z))
H0: \(\lambda_{ik}^{XZ}\neq0\) (Ada interaksi antara jenis kelamin
Tingkat Signifikansi
\(\alpha=5\%\)
Statistik Uji
\(\Delta\)Deviance = Deviance model conditional on Y - Deviance model homogenous = 3.3597 - 2.9849 = 0.3738
\(db\) = \(db\) model homogenous - \(db\) model saturated = 4 - 2 = 2
Daerah Penolakan
Tolak H0 jika \(\Delta\)Deviance \(\>>\) \(\chi_{0.05,db}^2 = \chi_{0.05,2}^2 = 5.991\)
Keputusan
Karena 0.3738 < 5.991, maka terima H0
Interpretasi
Dengan taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi antara jenis kelamin dan tingkat stress. Dengan kata lain, model yang terbentuk adalah model tanpa parameter \(\lambda_{ik}^{XZ}\)
# Deviance of Model
Deviance.model <- model_conditional_Y$deviance - model_homogenous$deviance # model_conditional_Y: conditional on Y, model_homogenous: homogenous
Deviance.model
## [1] 0.3737977
# Chi Square tabel dengan alpha = 0.05
derajat.bebas <- (4 - 2)
derajat.bebas
## [1] 2
chi.tabel <- qchisq((1 - 0.05), df = derajat.bebas)
chi.tabel
## [1] 5.991465
Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima", "Tolak")
Keputusan
## [1] "Terima"
Interpretasi Karena nilai Deviance.model = 0.3738 lebih kecil dari nilai kritis chi-square tabel = 5.99 (dengan df = 2, alpha = 0.05), maka keputusan uji adalah “Terima”.
Pada taraf nyata 5%, belum cukup bukti untuk menolak H0, atau dengan kata lain tidak ada interaksi antara jenis kelamin (X) dan tingkat stress (Z). Model yang terbentuk cukup hanya sampai dua interaksi dengan Y (conditional on Y), sehingga interaksi X*Z tidak signifikan secara statistik.
Model log-linear conditional pada Z memasukkan efek utama dan interaksi dua arah antara X dengan Z dan Y dengan Z, tanpa interaksi antara X dengan Y maupun interaksi tiga arah. \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} \]
# Conditional Association on Z
model_conditional_Z <- glm(counts ~ x.sex + y.freq + z.stres +
x.sex*z.stres + y.freq*z.stres,
family = poisson(link = "log"))
summary(model_conditional_Z)
##
## Call:
## glm(formula = counts ~ x.sex + y.freq + z.stres + x.sex * z.stres +
## y.freq * z.stres, family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.3401 0.1092 39.730 < 2e-16 ***
## x.sex1L -0.1613 0.1521 -1.060 0.289
## y.freq1jarang -1.4901 0.1957 -7.615 2.64e-14 ***
## z.stres1tinggi -1.4230 0.2147 -6.627 3.43e-11 ***
## z.stres2sedang -0.7010 0.1794 -3.907 9.35e-05 ***
## x.sex1L:z.stres1tinggi 0.1072 0.2240 0.479 0.632
## x.sex1L:z.stres2sedang 0.1041 0.2274 0.458 0.647
## y.freq1jarang:z.stres1tinggi 2.6251 0.2739 9.585 < 2e-16 ***
## y.freq1jarang:z.stres2sedang 1.3757 0.2588 5.316 1.06e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 124.9856 on 11 degrees of freedom
## Residual deviance: 3.0566 on 3 degrees of freedom
## AIC: 85.239
##
## Number of Fisher Scoring iterations: 4
Pengujian Ada Tidaknya Interaksi Antara X dan Y (Homogenous Model vs Conditional Association on Z)
Hipotesis
H0: \(\lambda_{ij}^{XY}=0\) (Tidak ada interaksi antara jenis kelamin (X) dan frekuensi olahraga (Y)) -
H0: \(\lambda_{ij}^{XY}\neq0\) (Ada interaksi antara jenis kelamin
Tingkat Signifikansi
\(\alpha=5\%\)
Statistik Uji
\(\Delta\)Deviance = Deviance model conditional on Y - Deviance model homogenous = 3.0566 - 2.9849 = 0.072
\(db\) = \(db\) model homogenous - \(db\) model saturated = 3 - 2 = 1
Daerah Penolakan
Tolak H0 jika \(\Delta\)Deviance \(\>>\) \(\chi_{0.05,db}^2 = \chi_{0.05,1}^2 = 3.841\)
Keputusan
Karena 0.072 < 3.841, maka tolak H0
Interpretasi
Dengan taraf nyata 5%, menolak H0 atau dapat dikatakan bahwa ada interaksi antara jenis kelamin dan frekuensi olahraga. Dengan kata lain, model yang terbentuk adalah model parameter \(\lambda_{ij}^{XY}\)
# Deviance of Model
Deviance.model <- model_conditional_Z$deviance - model_homogenous$deviance # model_conditional_Z: conditional on Z, model_homogenous: homogenous
Deviance.model
## [1] 0.07173518
# Chi Square tabel dengan alpha = 0.05
derajat.bebas <- (3 - 2)
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", "Tolak")
Keputusan
## [1] "Terima"
Interpretasi Karena nilai Deviance.model = 0.072 lebih kecil dari nilai kritis chi-square tabel = 3.84 (dengan df = 1, alpha = 0.05), maka keputusan uji adalah “Terima H0”.
Pada taraf nyata 5%, belum cukup bukti untuk menolak H0, atau dengan kata lain tidak ada interaksi antara jenis kelamin (X) dan frekuesni olahraga (Y). Model yang terbentuk cukup hanya sampai dua interaksi dengan Z (conditional on Z), sehingga interaksi XY tidak signifikan secara statistik.
| Model | Parameter | Deviance | Jumlah Parameter | df | AIC |
|---|---|---|---|---|---|
| Saturated | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ + λᵢⱼₖˣʸᶻ | 0.000 | 12 | 0 | 88.18 |
| Homogenous | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ + λⱼₖʸᶻ | 2.9849 | 10 | 2 | 87.168 |
| Conditional on X | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λᵢₖˣᶻ | 115.77 | 8 | 4 | 195.95 |
| Conditional on Y | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢⱼˣʸ + λⱼₖʸᶻ | 3.3597 | 8 | 4 | 83.541 |
| Conditional on Z | λ + λᵢˣ + λⱼʸ + λₖᶻ + λᵢₖˣᶻ + λⱼₖʸᶻ | 3.0566 | 9 | 3 | 85.239 |
| Interaksi | Pengujian | \(\Delta\)Deviance | \(\Delta\)df | Chi-Square Tabel | Keputusan | Keterangan |
|---|---|---|---|---|---|---|
| XYZ | Saturated vs Homogenous | 2.9849 | 2 | 5.991 | Terima H0 | tidak ada interaksi |
| YZ | Conditional on X vs Homogenous | 112.7855 | 2 | 5.991 | Tolak H0 | ada interaksi |
| XZ | Conditional on Y vs Homogenous | 0.3738 | 2 | 5.991 | Terima H0 | tidak ada interaksi |
| XY | Conditional on Z vs Homogenous | 0.072 | 1 | 3.841 | Terima H0 | tidak ada interaksi |
Dari hasil di atas diketahui bahwa asosiasi yang nyata hanya terdapat antara frekuensi olahraga dan tingkat stress. Sehingga, model terbaik adalah:
\[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{YZ}_{jk} \]
Model terbaik adalah model log-linear tanpa interaksi tiga arah dan hanya memuat interaksi dua arah antara frekuesni olahraga (Y) dan tingkat stress (Z).
Model terbaik dipilih berdasarkan pengujian interaksi yang signifikan, yaitu hanya interaksi dua arah antara frekuensi olahraga (Y) dan tingkat stress (Z):
\[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{YZ}_{jk} \]
# Model Terbaik
bestmodel <- glm(counts ~ x.sex + y.freq + z.stres +
y.freq*z.stres,
family = poisson(link = "log"))
summary(bestmodel)
##
## Call:
## glm(formula = counts ~ x.sex + y.freq + z.stres + y.freq * z.stres,
## family = poisson(link = "log"))
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 4.30920 0.09492 45.398 < 2e-16 ***
## x.sex1L -0.09531 0.09315 -1.023 0.306
## y.freq1jarang -1.49009 0.19568 -7.615 2.64e-14 ***
## z.stres1tinggi -1.37231 0.18660 -7.354 1.92e-13 ***
## z.stres2sedang -0.65176 0.14337 -4.546 5.47e-06 ***
## y.freq1jarang:z.stres1tinggi 2.62507 0.27386 9.585 < 2e-16 ***
## y.freq1jarang:z.stres2sedang 1.37568 0.25876 5.316 1.06e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 124.9856 on 11 degrees of freedom
## Residual deviance: 3.3587 on 5 degrees of freedom
## AIC: 81.541
##
## Number of Fisher Scoring iterations: 4
| Model | AIC |
|---|---|
| Saturated | 88.18 |
| Homogenous | 87.168 |
| Conditional on X | 195.95 |
| Conditional on Y | 83.541 |
| Conditional on Z | 85.239 |
| Model Terbaik | 81.541 |
Dari hasil di atas, terlihat bahwa model terbaik memiliki AIC yang lebih rendah dibandingkan dengan model lainnya seperti saturated, homogennous, dan conditional model.
# Interpretasi koefisien model terbaik
data.frame(
koef = bestmodel$coefficients,
exp_koef = exp(bestmodel$coefficients)
)
| Koefisien | Nilai Koefisien | Exp(Nilai Koefisien) |
|---|---|---|
| Intercept | 4.309 | 74.3809 |
| x.sex1M | -0.0953 | 0.9091 |
| y.freq1jarang | -1.49 | 0.225 |
| z.stress1tinggi | -1.372 | 0.2535 |
| z.stress2sedang | -0.652 | 0.5211 |
| z.stress1tinggi:y.freq1jarang | 2.625 | 13.806 |
\(\exp(\lambda_{\text{1M}}^X) = \exp(-0.0953) = 0.909 \rightarrow\) nilai odds
Tanpa memperhatikan frekuensi olahraga dan tingkat stres, peluang
seseorang berjenis kelamin laki-laki adalah
0,909 kali dibandingkan
perempuan.
Atau, peluang seseorang perempuan adalah \(1 / 0.909 = 1.10\) kali dibandingkan
laki-laki.
\(\exp(\lambda_{\text{1jarang}}^Y) = \exp(-1.49) = 0.225 \rightarrow\) nilai odds
Tanpa memperhatikan jenis kelamin dan tingkat stres, peluang
seseorang yang berolahraga jarang adalah 0.225
kali dibandingkan yang sering.
Atau, peluang seseorang yang sering berolahraga adalah
\(1 / 0.225 = 4.44\) kali dibandingkan
yang jarang.
\(\exp(\lambda_{\text{1tinggi}}^Z) = \exp(-1.372) = 0.253 \rightarrow\) nilai odds
Tanpa memperhatikan jenis kelamin dan frekuensi olahraga, peluang
seseorang dengan tingkat stres tinggi adalah
0.253 kali dibandingkan yang tingkat stresnya
rendah.
Atau, peluang seseorang dengan stres rendah adalah
\(1 / 0.253 = 3.95\) kali dibandingkan
dengan yang stres tinggi.
\(\exp(\lambda_{\text{2sedang}}^Z) = \exp(-0.652) = 0.521 \rightarrow\) nilai odds
Tanpa memperhatikan jenis kelamin dan frekuensi olahraga, peluang
seseorang dengan tingkat stres sedang adalah
0.521 kali dibandingkan yang
rendah.
Atau, peluang seseorang dengan stres rendah adalah
\(1 / 0.521 = 1.92\) kali dibandingkan
yang stres sedang.
\(\exp(\lambda_{\text{1tinggi,1jarang}}^{Z,Y}) = \exp(2.625) = 13.806 \rightarrow\) nilai odds ratio
Tanpa memperhatikan jenis kelamin, odds seseorang yang memiliki stres tinggi dan olahraga jarang adalah 13.806 kali lebih tinggi dibandingkan dengan orang yang tidak memiliki kedua kondisi tersebut.
Secara manual, nilai fitted value diperoleh dengan cara sebagai berikut:
\[\begin{align*} \\ \hat{\mu}_{111} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{1jarang} + \lambda^{z}_{1tinggi} + \lambda^{yz}_{1jarang,1tinggi} \right) \\ \\ &= \exp\left( 4.309 - 0.0953 - 1.49 - 1.372 + 2.625 \right) \\ \\ &= \exp(3.9767) = 53.333 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{112} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{1jarang} + \lambda^{z}_{2sedang} + \lambda^{yz}_{1jarang,2sedang} \right) \\ \\ &= \exp\left( 4.309 - 0.0953 - 1.49 - 0.652 + 1.376 \right) \\ \\ &= \exp(3.4477) = 31.429 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{113} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{1jarang} + \lambda^{z}_{3rendah} + \lambda^{yz}_{1jarang,3rendah} \right) \\ \\ &= \exp\left( 4.309 - 0.0953 - 1.49 + 0 + 0 \right) \\ \\ &= \exp(2.7237) = 15.2381 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{121} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{2rutin} + \lambda^{z}_{1tinggi} + \lambda^{yz}_{2rutin,1tinggi} \right) \\ \\ &= \exp\left( 4.309 - 0.0953 + 0 - 1.372 + 0 \right) \\ \\ &= \exp(2.8417) = 17.143 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{122} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{2rutin} + \lambda^{z}_{2sedang} + \lambda^{yz}_{2rutin,2sedang} \right) \\ \\ &= \exp\left( 4.309 - 0.0953 + 0 - 0.652 + 0 \right) \\ \\ &= \exp(3.5617) = 35.238 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{123} &= \exp\left( \lambda + \lambda^{x}_{1m} + \lambda^{y}_{2rutin} + \lambda^{z}_{3rendah} + \lambda^{yz}_{2rutin,3rendah} \right) \\ \\ &= \exp\left( 4.309 - 0.0953 + 0 + 0 + 0 \right) \\ \\ &= \exp(4.2137) = 67.619 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{211} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{1jarang} + \lambda^{z}_{1tinggi} + \lambda^{yz}_{1jarang,1tinggi} \right) \\ \\ &= \exp\left( 4.309 + 0 - 1.49 - 1.372 + 2.625 \right) \\ \\ &= \exp(4.072) = 58.667 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{212} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{1jarang} + \lambda^{z}_{2sedang} + \lambda^{yz}_{1jarang,2sedang} \right) \\ \\ &= \exp\left( 4.309 + 0 - 1.49 - 0.652 + 1.376 \right) \\ \\ &= \exp(3.543) = 34.571 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{213} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{1jarang} + \lambda^{z}_{3rendah} + \lambda^{yz}_{1jarang,3rendah} \right) \\ \\ &= \exp\left( 4.309 + 0 - 1.49 + 0 + 0 \right) \\ \\ &= \exp(2.819) = 16.762 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{221} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{2rutin} + \lambda^{z}_{1tinggi} + \lambda^{yz}_{2rutin,1tinggi} \right) \\ \\ &= \exp\left( 4.309 + 0 + 0 - 1.372 + 0 \right) \\ \\ &= \exp(2.937) = 18.857 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{222} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{2rutin} + \lambda^{z}_{2sedang} + \lambda^{yz}_{2rutin,2sedang} \right) \\ \\ &= \exp\left( 4.309 + 0 + 0 - 0.652 + 0 \right) \\ \\ &= \exp(3.657) = 38.762 \\ \end{align*}\]
\[\begin{align*} \\ \hat{\mu}_{223} &= \exp\left( \lambda + \lambda^{x}_{2f} + \lambda^{y}_{2rutin} + \lambda^{z}_{3rendah} + \lambda^{yz}_{2rutin,3rendah} \right) \\ \\ &= \exp\left( 4.309 + 0 + 0 + 0 + 0 \right) \\ \\ &= \exp(4.309) = 74.381 \\ \end{align*}\]
Keterangan: Nilai \(\hat{\mu}_{ijk}\) akan sama apapun referensi dari kategori peubahnya yang digunakan.
Dengan menggunakan R:
# Fitted values dari model terbaik
data.frame(
Stress = z.stres,
sex = x.sex,
Freq = y.freq,
counts = counts,
fitted = bestmodel$fitted.values
)
BBC Indonesia. (2018, Juli 11). Ilmuwan sebut perempuan lebih merasakan efek alkohol daripada laki-laki. BBC News Indonesia. Diakses dari https://www.bbc.com/indonesia/vert-fut-44789893.
Damayanti, A., & Umeda, M. (2019). Hubungan konsumsi alkohol dengan perilaku seks bebas pada remaja. Indonesian Journal of Nursing Sciences and Practices, 2(1), 21-26
Modul Analasis Data Kategori Pertemuan 14
Regina Septie Nugraheni. (n.d.). BAB IV. Universitas Katolik Soegijapranata. Retrieved from https://repository.unika.ac.id/15599/5/13.70.0074%20Regina%20Septie%20Nugraheni%20BAB%20IV.pdf
Rizka Putri. (n.d.). Body dissatisfaction dan eating behavior pada remaja putri. Universitas Islam Negeri Syarif Hidayatullah Jakarta. Retrieved from https://repository.uinjkt.ac.id/dspace/bitstream/123456789/64420/1/11181330000022_RIZKA%20PUTRI%20-%20Rizka%20Putri.pdf
Unknown Author. (2020). Pengaruh body dissatisfaction terhadap perilaku diet remaja putri. Seminar Nasional Pendidikan Nusantara, STKIP Kusumanegara. Retrieved from https://jurnal.stkipkusumanegara.ac.id/index.php/semnara2020/article/download/844/349
Unknown Author. (n.d.). Hubungan antara body dissactisfaction dengan kecenderungan perilaku diet pada remaja putri. Garuda.kemdikbud.go.id. Retrieved from http://download.garuda.kemdikbud.go.id/article.php?article=1066894&val=8723&title=Hubungan%20Antara%20Body%20Dissactisfaction%20Dengan%20Kecenderungan%20Perilaku%20Diet%20Pada%20Remaja%20Putri
Unknown Author. (n.d.). Artikel penelitian di CORISINDO (UTB). Retrieved from https://corisindo.utb-univ.ac.id/index.php/penelitian/article/download/120/56/1119
Unknown Author. (n.d.). OSF Preprint [PDF]. Retrieved from https://osf.io/v9a7x/download/?format=pdf
Unknown Author. (n.d.). JTIKOR – Jurnal Teknik Informatika dan Komputer. Universitas Pendidikan Indonesia. Retrieved from https://ejournal.upi.edu/index.php/JTIKOR/article/view/24895