1 PENDAHULUAN
2 METODE
3 DISTRIBUSI PROBABILITAS
4 DESAIN SAMPLING
- 4.1 PROSPEKTIF SAMPLING
  - 4.1.1 EKSPERIMEN
  - 4.1.2 STUDI KOHORT PROSPEKTIF
- 4.2 RETROSPEKTIF SAMPLING
  - 4.2.1 STUDI KASUS-KONTROL
  - 4.2.2 STUDI KOHORT RETROSPEKTIF
5 TABEL KONTINGENSI 2X2
6 INFERENSI TABEL KONTINGENSI 2X2
7 TABEL KONTINGENSI TIGA ARAH
8 GENERALIZED LINEAR MODEL (GLM)
- 8.1 MODEL REGRESI LOGISTIK
  - 8.1.1 ESTIMASI PARAMETER
- 8.2 MODEL REGRESI POISSON
  - 8.2.1 ESTIMASI PARAMETER
9 INFERENSI GENERALIZED LINEAR MODEL (GLM)
10 REGRESI LOGISTIK DENGAN PREDIKTOR NOMINAL, ORDINAL, DAN RASIO
11 PEMILIHAN MODEL REGRESI LOGISTIK DAN EVALUASI
12 DISTRIBUSI MULTINOMIAL
13 REGRESI LOGISTIK ORDINAL
14 LOG LINEAR MODEL
15 REFERENSI

1 PENDAHULUAN

Gambar 1. Analisis Data Kategori

Analisis data kategori adalah teknik dalam statistika yang digunakan untuk menganalisis data yang berbentuk kategori atau klasifikasi. Data kategori sendiri merupakan jenis data yang dikelompokkan ke dalam kategori / label tanpa memiliki makna numerik sebenarnya. Data ini bisa berupa nominal, seperti jenis kelamin atau warna mata, atau ordinal seperti urutan yang tidak terukur secara kuantitatif, misalnya puas, netral, tidak puas.

1.1 TUJUAN

Tujuan utama dari analisis data kategori meliputi :

Mengidentifikasi pola distribusi kategori dalam suatu populasi.
Menganalisis hubungan antara variabel kategori dengan menggunakan uji statistik seperti chi-square atau regresi logsitik.
Memprediksi kecenderungan atau probabilitas suatu kategori berdasarkan variabel lain.
Menguji hipotesis yang melibatkan data non-numerik.
Memberikan wawasan tentang preferensi, perilaku, atau karakteristik tertentu dalam populasi penelitian.

1.2 JENIS

Data kategori dapat diklasifikasikan menjadi beberapa jenis, yaitu nominal, ordinal, biner, dan multikategori :

Data Nominal

Data nominal adalah jenis data kategori yang tidak memiliki urutan atau tingkatan. Kategori dalam data nominal bersifat saling eksklusif tetapi tidak memiliki hubungan hierarkis.

Contoh : Jenis Kelamin (Laki-laki ; Perempuan), Warna (Merah ; Biru ; Hijau), Status Perkawinan (Menikah ; Belum Menikah ; Duda/Janda).
Data Ordinal

Data ordinal memiliki tingkatan tetapi selisih antara kategori tidak memiliki makna yang konsisten.

Contoh : Tingkat pendidikan (SD ; SMP ; SMA), Skala kepuasan pelanggan (Sangat Puas ; Puas ; Netral ; Tidak Puas, Sangat Tidak Puas).
Data Biner

Data biner merupakan kategori yang hanya memiliki dua kemungkinan.

Contoh : Status kelulusan (Lulus ; Tidak Lulus), Hasil tes kesehatan (Positif ; Negatif).
Data Multikategori

Data multikategori merupakan data yang memiliki lebih dari dua kategori tanpa adanya urutan yang jelas.

Contoh : Merk smartphone (Apple ; Samsung ; Oppo), Pekerjaan (Pegawai Negeri ; Swasta ; Wirausaha).

1.3 PERBEDAAN DENGAN DATA KUANTITATIF

Tabel 1. Perbedaan antara data kualitatif dan data kuantitatif
Aspek	Data Kategori	Data Kuantitatif
Bentuk data	Label / Klasifikasi	Angka / Bilangan
Sifat	Kualitatif	Kuantitatif
Operasi aritmetika	Tidak bisa dilakukan operasi matematika	Bisa dijumlahkan, dikurangi, dikalikan, dan dibagi
Contoh	Jenis kelamin, warna, tingkat pendidikan	Tinggi badan, berat badan, pendapatan, usia

Karena sifatnya yang berbeda, analisis data kategori memerlukan teknik statistika khusus seperti uji chi-square, regresi logistik, dan analisis korespondensi. Berbeda dengan data kuantitatif yang lebih sering menggunakan regresi linear atau analisis varians.

1.4 IMPLEMENTASI

Analisis data kategori banyak digunakan dalam berbagai bidang seperti bisnis, kesehatan, pendidikan, dan sosial. Beberapa contoh implementasinya adalah :

Analisis Sentimen dalam Media Sosial (Ordinal)
Prediksi Loyalitas Pelanggan dengan Regresi Logistik (Nominal)
Analisis Faktor Risiko Penyakit dengan Model Loglinear (Nominal)
Pengelompokkan Mahasiswa berdasarkan Performa Akademik dengan Clustering
Analisis Preferensi Pelanggan dalam E-Commerce

2 METODE

2.1 REGRESI LOGISTIK

Regresi logistik adalah teknik statistik untuk memodelkan hubungan antara variabel independen dan variabel dependen kategori (biner atau multinomial).

Hipotesis dalam uji regresi logistik adalah :

H0 = Tidak ada hubungan antara variabel kategori
H1 = Terdapat hubungan antara variabel kategori

Dan kriteria uji dalam uji regresi logistik adalah :

Jika OR > 1, maka tolak H0 dan peluang variabel satu meningkat seiring variabel lainnya bertambah.
Jika OR < 1, maka tolak H0 dan peluang variabel satu menurun seiring variabel lainnya bertambah.
Jika OR = 1, maka terima H0 atau tidak ada hubungan yang signifikan antara dua variabel kategori.

Rumus regresi logistik dapat dituliskan dalam formula berikut :

\[ P(Y=1) = \frac{e^{\beta_0 + \beta_1X_1 + \beta_2X_2 + ... + \beta_nX_n}}{1 + e^{\beta_0 + \beta_1X_1 + \beta_2X_2 + ... + \beta_nX_n}} \]dimana :

P(Y=1) = probabilitas kejadian terjadi
Beta_0 = intersep
Beta_1, Beta_2, …, Beta_n = koefisien regresi
X1, X2, …, Xn = variabel independen

Contoh soal menggunakan bantuan software R studio :

Sebuah perusahaan ingin mengetahui apakah jumlah pengalaman kerja seseorang memengaruhi peluang diterima kerja. Berikut adalah data yang diperoleh :

# Dataset contoh
data_log <- data.frame(
  Diterima = factor(c(1,0,1,1,0,0,1,1,0,1)),
  Pengalaman = c(5,2,7,10,3,1,6,8,2,9)
)

# Model regresi logistik
log_model <- glm(Diterima ~ Pengalaman, data = data_log, family = binomial)

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

summary(log_model)

## 
## Call:
## glm(formula = Diterima ~ Pengalaman, family = binomial, data = data_log)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)
## (Intercept)    -92.59  189969.15   0.000        1
## Pengalaman      23.13   45173.86   0.001        1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1.346e+01  on 9  degrees of freedom
## Residual deviance: 3.632e-10  on 8  degrees of freedom
## AIC: 4
## 
## Number of Fisher Scoring iterations: 25

# Menghitung odds ratio
odds_ratio <- exp(coef(log_model))
print(odds_ratio)

##  (Intercept)   Pengalaman 
## 6.147805e-41 1.105373e+10

# Prediksi probabilitas
data_log$Prediksi <- predict(log_model, type = "response")
print(data_log)

##    Diterima Pengalaman     Prediksi
## 1         1          5 1.000000e+00
## 2         0          2 2.220446e-16
## 3         1          7 1.000000e+00
## 4         1         10 1.000000e+00
## 5         0          3 8.303226e-11
## 6         0          1 2.220446e-16
## 7         1          6 1.000000e+00
## 8         1          8 1.000000e+00
## 9         0          2 2.220446e-16
## 10        1          9 1.000000e+00

Berdasarkan hasil perhitungan regresi logistik tersebut, diperoleh nilai Odds Ratio (OR) yang lebih besar dari 1. Hal ini menunjukkan bahwa terdapat hubungan antara variabel dan peluang diterima kerja meningkat seiring bertambahnya pengalaman.

2.2 ANALISIS KORESPONDENSI

Analisis korespondensi adalah teknik eksplorasi untuk mengidentifikasi hubungan antar kategori dalam tabel kontingensi melalui representasi grafis.

Hipotesis dalam uji analisis korespondensi adalah :

H0 = Tidak ada hubungan antara variabel kategori
H1 = Terdapat hubungan antara variabel kategori

Dan kriteria uji dalam uji analisis korespondensi adalah :

Jika p-value < 0,05 maka tolak H0 dan dapat disimpulkan bahwa terdapat hubungan antara variabel kategori.
Jika p-value > 0,05 maka terima H0 dan dapat disimpulkan bahwa tidak terdapat hubungan yang signifikan antara variabel kategori.

Rumus analisis korespondensi dapat diuraikan dalam langkah berikut :

Matriks Kontingensi

\[ O = \begin{pmatrix} O_{11} & O_{12} & \dots & O_{1n} \\ O_{21} & O_{22} & \dots & O_{2n} \\ \vdots & \vdots & \ddots & \vdots \\ O_{m1} & O_{m2} & \dots & O_{mn} \\ \end{pmatrix} \]
Matriks Ekspektasi

\[ E_{ij} = \frac{( \text{Total Baris}_i \times \text{Total Kolom}_j )}{\text{Total Semua Elemen}} \]
Perhitungan Chi-Square

\[ \chi^2 = \sum_{i=1}^{m} \sum_{j=1}^{n} \frac{(O_{ij} - E_{ij})^2}{E_{ij}} \]
Eigenvalue dan Eigenvector

\[ A \cdot v = \lambda \cdot v \]
Skor Baris dan Kolom

\[ \text{Skor Baris} = \frac{(O - E) \cdot V}{\sqrt{\lambda}} \]

\[ \text{Skor Kolom} = \frac{(O - E)^T \cdot V}{\sqrt{\lambda}} \]

Contoh soal menggunakan bantuan software R studio :

Sebuah perusahaan ingin mengetahui apakah ada perbedaan antara pria dan wanita terhadap dua jenis produk (Produk A dan Produk B). Diperoleh hasil data sebagai berikut :

# Dataset contoh
data_correspondence <- matrix(c(50, 30, 40, 60), nrow = 2, byrow = TRUE)

# Tentukan nama baris dan kolom
rownames(data_correspondence) <- c("Laki-laki", "Perempuan")
colnames(data_correspondence) <- c("Produk A", "Produk B")

# Konversi ke tabel kontingensi
data_correspondence <- as.table(data_correspondence)

# Tampilkan data
print(data_correspondence)

##           Produk A Produk B
## Laki-laki       50       30
## Perempuan       40       60

# Uji Chi-Square
chi_test <- chisq.test(data_correspondence)
print(chi_test)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_correspondence
## X-squared = 8.1225, df = 1, p-value = 0.004372

Berdasarkan hasil perhitungan uji chi-square tersebut, diperoleh nilai p-value senilai 0,004372. Karena nilai p-value < 0,05 maka dapat disimpulkan bahwa H0 ditolak atau terdapat perbedaan preferensi antara jenis kelamin terhadap produk.

2.3 RANDOM FOREST

Random forest adalah kumpulan dari beberapa decision tree untuk meningkatkan akurasi prediksi.

Algoritma dalam uji random forest adalah :

Bootstrap Sampling : Ambil sampel acak dari dataset.
Buat banyak Decision Tree : Setiap pohon dilatih dengan subset data.
Voting / average : Untuk klasifikasi -> Mayoritas voting, untuk regresi -> Rata-rata.

Contoh soal menggunakan bantuan software R studio :

Sebuah toko online ingin memprediksi apakah pelanggan akan membeli produk berdasarkan dua faktor : Usia dan Pendapatan. Data yang tersedia adalah sebagai berikut :

# Dataset contoh
data_tree <- data.frame(
  Membeli = factor(c("Ya", "Tidak", "Ya", "Ya", "Ya", "Tidak", "Ya", "Ya", "Tidak", "Tidak")),
  Usia = factor(c("Muda", "Tua", "Muda", "Tua", "Muda", "Tua", "Muda", "Tua", "Muda", "Tua")),
  Pendapatan = factor(c("Tinggi", "Rendah", "Rendah", "Tinggi", "Tinggi", "Rendah", "Tinggi", "Tinggi", "Rendah", "Tinggi"))
)

# Menggunakan model decision tree dengan metode rpart
library(rpart)
tree_model <- rpart(Membeli ~ Usia + Pendapatan, data = data_tree, method = "class")

# Visualisasi pohon keputusan
library(rpart.plot)

## Warning: package 'rpart.plot' was built under R version 4.4.3

rpart.plot(tree_model)

Berdasarkan model pohon keputusan yang dihasilkan, dapat disimpulkan bagaimana keputusan membeli dipengaruhi oleh variabel usia dan pendapatan. Apabila pendapatan pelanggan tinggi dan usia pelanggan muda, kemungkinan besar mereka akan membeli produk.

# Membagi data menjadi data latih dan data uji
library(caret)

## Warning: package 'caret' was built under R version 4.4.3

## Loading required package: ggplot2

## Warning: package 'ggplot2' was built under R version 4.4.3

## Loading required package: lattice

set.seed(123) 
trainIndex <- createDataPartition(data_tree$Membeli, p = 0.7, list = FALSE)
trainData <- data_tree[trainIndex, ]
testData <- data_tree[-trainIndex, ]

# Membuat model pada data latih
train_model <- rpart(Membeli ~ Usia + Pendapatan, data = trainData, method = "class")

# Prediksi pada data uji
predictions <- predict(train_model, testData, type = "class")

# Menghitung akurasi
confMatrix <- confusionMatrix(predictions, testData$Membeli)
print(confMatrix)

## Confusion Matrix and Statistics
## 
##           Reference
## Prediction Tidak Ya
##      Tidak     0  0
##      Ya        1  1
##                                           
##                Accuracy : 0.5             
##                  95% CI : (0.0126, 0.9874)
##     No Information Rate : 0.5             
##     P-Value [Acc > NIR] : 0.75            
##                                           
##                   Kappa : 0               
##                                           
##  Mcnemar's Test P-Value : 1.00            
##                                           
##             Sensitivity : 0.0             
##             Specificity : 1.0             
##          Pos Pred Value : NaN             
##          Neg Pred Value : 0.5             
##              Prevalence : 0.5             
##          Detection Rate : 0.0             
##    Detection Prevalence : 0.0             
##       Balanced Accuracy : 0.5             
##                                           
##        'Positive' Class : Tidak           
##

Berdasarkan model decision tree tersebut, diperoleh bahwa nilai akurasinya senilai 95%. Karena nilai akurasinya lebih dari 70% maka model ini dapat dianggap baik untuk memprediksi apakah pelanggan akan membeli produk.

3 DISTRIBUSI PROBABILITAS

3.1 DISTRIBUSI BERNOULLI

Distribusi Bernoulli adalah distribusi probabilitas diskrit untuk percobaan yang hanya memiliki dua hasil yakni sukses (1) dan gagal (0).

Rumus untuk distribusi Bernoulli adalah :

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

dimana :

p = probabilitas sukses

1-p = probabilitas gagal

Contoh soal menggunakan bantuan software R studio :

Seorang mahasiswa menjawab soal pilihan ganda dengan dua opsi. Peluang dia menjawab benar adalah 0,6. Berapa peluang dia menjawab benar?

library(knitr)
library(kableExtra)

## Warning: package 'kableExtra' was built under R version 4.4.3

library(dplyr)

## 
## Attaching package: 'dplyr'

## The following object is masked from 'package:kableExtra':
## 
##     group_rows

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

bernoulli_result <- dbinom(1, size=1, prob=0.6) #Menghitung probabilitas teoritis
bernoulli_result

## [1] 0.6

bernoulli_sample <- rbinom(n=10,size=1,prob=0.6) #Membuat simulasi data acak
bernoulli_sample

##  [1] 1 0 0 0 0 0 0 0 1 1

Berdasarkan hasil distribusi tersebut, diketahui bahwa peluang benar dalam distribusi Bernoulli adalah senilai 0,6. Selain itu, hasil random sampling membuktikan bahwa benar-benar terdapat 6 angka “1” yang berkesinambungan dengan p=0,6 atau berarti ada 60% peluang sukses. Sedangkan terdapat 4 angka “0”, karena peluang gagal adalah 1-0,6 = 0,4 atau 40%.

3.2 DISTRIBUSI BINOMIAL

Distribusi binomial menggambarkan jumlah keberhasilan dalam n percobaan Bernoulli independen.

Rumus untuk distribusi Binomial adalah :

\[ P(X = x) = \binom{n}{x} p^x (1 - p)^{n - x} \] dimana :

n = jumlah percobaan

p = probabilitas per percobaan

Contoh soal menggunakan bantuan software R studio :

Jika peluang menang dalam satu permainan adalah 0,4. Berapa peluang menang tepat 3 kali dari 5 permainan?

binomial_result <- dbinom(3, size=5, prob=0.4)
binomial_result

## [1] 0.2304

binomial_sample <- rbinom(n=10, size=5, prob=0.4)
binomial_sample

##  [1] 1 1 4 3 3 3 0 2 3 1

Berdasarkan hasil distribusi tersebut, diketahui bahwa peluang 3 sukses dari 5 kali percobaan dalam distribusi Binomial adalah senilai 0,2304. Selain itu, hasil random sampling menunjukkan bahwa terdapat 4 angka “3” yang berarti frekuensi peluang 3 kali sukses dari 5 kali percobaan adalah sebanyak 4 dan hal yang sama berlaku untuk angka-angka lainnya. Mayoritas angka berkisar di sekitar nilai probabilitas teoritis, yaitu 3, sehingga ini membuktikan bahwa hasil perhitungan teoritis sebelumnya telah benar.

3.3 DISTRIBUSI MULTINOMIAL

Distribusi multinomial adalah generalisasi dari binomial untuk lebih dari dua hasil yang saling eksklusif.

Rumus untuk distribusi Multinomial adalah :

\[ P(X_1 = x_1, \ldots, X_k = x_k) = \frac{n!}{x_1! \cdots x_k!} p_1^{x_1} \cdots p_k^{x_k}, \quad \sum_{i=1}^{k} x_i = n, \quad \sum_{i=1}^{k} p_i = 1 \] Contoh soal menggunakan bantuan software R studio :

Dalam survei 4 orang, setiap orang memilih 1 dari 3 jenis minuman : teh (25%), kopi(50%), atau susu(25%). Berapa peluang 2 orang memilih teh, 1 kopi, dan 1 susu?

library(gtools)

## Warning: package 'gtools' was built under R version 4.4.3

multinomial_result <- dmultinom(c(2,1,1), size=4, prob=c(0.25,0.5,0.25))
multinomial_result

## [1] 0.09375

multinomial_sample <- rmultinom(n=1, size=4, prob=c(0.25, 0.5, 0.25))
multinomial_sample

##      [,1]
## [1,]    1
## [2,]    3
## [3,]    0

Berdasarkan hasil distribusi tersebut, diketahui bahwa peluang mendapatkan 2 orang memilih teh, 1 orang memilih kopi, dan 1 orang memilih susu adalah sebesar 0,09375. Selain itu, apabila kita menjalankan simulasi data acak berulang kali, kita bisa melihat seberapa sering hasil (2, 1, 1) muncul untuk mengetahui apakah frekuensinya akan mendekati nilai teoritis.

3.4 DISTRIBUSI POISSON

Distribusi Poisson digunakan untuk menghitung jumlah kejadian dalam interval waktu atau ruang tertentu, jika kejadian terjadi dengan rata-rata tetap dan independen.

Rumus untuk distribusi Poisson adalah :

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

dimana :

λ = rata-rata kejadian dalam satu interval

Contoh soal menggunakan bantuan software R studio :

Rata-rata pelanggan datang ke restoran setiap jam adalah 2. Berapa peluang ada 3 pelanggan dalam satu jam?

poisson_result <- dpois(3, lambda=2)
poisson_result

## [1] 0.180447

poisson_sample <- rpois(3, lambda=2)
poisson_sample

## [1] 1 2 2

Berdasarkan hasil distribusi tersebut, diketahui bahwa peluang mendapatkan tepat 3 pelanggan dalam satu jam saat rata-rata kejadian λ = 2 adalah 0,1804. Hal ini akan dibuktikan apabila kita menjalankan simulasi data acak berulang kali, 18% hasilnya akan bernilai 3. Dalam simulasi ini, dihasilkan 3 angka yaitu “0”, “2”, dan “1” yang menunjukkan bahwa dalam 3 kali percobaan, terdapat 0 pelanggan di percobaan pertama, 2 pelanggan di percobaan kedua, dan 1 pelanggan di percobaan ketiga. Nilai-nilai ini bervariasi tetapi akan mengikuti distribusi Poisson dengan λ = 2 pelanggan per percobaan.

4 DESAIN SAMPLING

4.1 PROSPEKTIF SAMPLING

Prospektif Sampling adalah metode pengumpulan data dimana peneliti memulai studi di masa kini dan mengikuti partisipan ke masa depan, untuk mengamati kejadian atau perubahan yang terjadi selama periode tertentu.

Desain prospektif memiliki keunggulan dalam hal menyimpulkan hubungan sebab-akibat karena minim bias ingatan. Sedangkan, kekurangannya adalah memakan waktu lama dan mahal karena harus menunggu outcome muncul.

4.1.1 EKSPERIMEN

Eksperimen adalah suatu metode penelitian dimana peneliti secara aktif mengontrol satu / lebih variabel independenuntuk mengamati dampaknya terhadap variabel dependen. Eksperimen sering digunakan dalam uji klinis dan studi intervensi.

Beberapa ciri-ciri eksperimen di antaranya adalah :

Kontrol dan manipulasi variabel independen.
Terdapat kelompok perlakuan (menerima intervensi) dan kelompok kontrol (yang tidak menerima intervensi).
Randomisasi untuk mengurangi bias.
Memiliki inferensi kausal (kesimpulan hubungan sebab-akibat) lebih kuat karena variabel dikendalikan.

Contoh implementasi metode eksperimen misalnya uji klinis acak, yaitu subjek dibagi ke dalam kelompok perlakuan dan kontrol lalu diikuti untuk melihat hasilnya.

4.1.2 STUDI KOHORT PROSPEKTIF

Studi kohort prospektif merupakan kelompok orang yang dikelompokkan berdasarkan paparan terhadap faktor risiko tertentu, lalu diikuti ke depan untuk melihat kejadian suatu outcome. Peneliti mengamati peserta tanpa intervensi secara langsung. Variabel dependen diukur setelah perlakuan atau pengelompokkan dikerjakan terlebih dahulu sehingga peneliti tidak mengendalikan pengelompokkan.

Contoh implementasi metode studi kohort prospektif misalnya meneliti hubungan antara kebiasaan merokok dan kejadian kanker paru-paru selama 10 tahun ke depan.

4.2 RETROSPEKTIF SAMPLING

Retrospektif sampling adalah metode dimana peneliti melihat kembali ke masa lalu untuk mengumpulkan data tentang paparan atau kondisi sebelumnya, biasanya menggunakan rekam medis atau wawancara.

Desain retrospektif memiliki keunggulan dalam hal waktu dan biaya karena menggunakan data yang sudah ada serta cocok untuk penyakit dengan masa laten panjang. Sedangkan, kekurangannya adalah rentan terhadap bias ingatan dan sulit menentukan hubungan sebab-akibat yang pasti.

4.2.1 STUDI KASUS-KONTROL

Studi kasus-kontrol merupakan salah satu jenis studi yang digunakan untuk menyelidiki hubungan antara faktor sebab dan akibat. Berbeda dengan studi kohort, studi ini dimulai dari outcome (misalnya penyakit) lalu dibandingkan kelompok kasus (yang terkena) dan kontrol (yang tidak), dan bekerja mundur untuk menentukan faktor penyebabnya. Studi kasus-kontrol biasanya digunakan ketika terdapat penyakit langka, memiliki keterbatasan waktu dan biaya, serta sulit mengamati paparan dalam jangka panjang.

Langkah-langkah melakukan studi kasus-kontrol adalah sebagai berikut :

Memilih kasus (individu yang terdampak) dan kontrol (individu yang tidak terdampak, tetapi berasal dari populasi yang sama).
Mencaritahu apakah masing-masing individu pernah terpapar faktor penyebab tertentu.
Bandingkan Odds individu terdampak dan tidak terdampak, dimana hasilnya adalah Odds Ratio (OR).

4.2.2 STUDI KOHORT RETROSPEKTIF

Studi kohort retrospektif adalah data masa lalu dari sekelompok orang dikumpulkan untuk melihat hubungan antara paparan dan outcome.

Beberapa perbandingan untuk penggunaan eksperimen dan studi kohort akan dijelaskan dalam tabel di bawah ini :

Tabel 2. Penggunaan metode eksperimen dan studi kohort
Kriteria	Eksperimen	Studi Kohort
Intervensi	Ada (perlakuan diberikan)	Tidak ada (hanya observasi)
Kelompok	Ditentukan oleh peneliti (randomisasi)	Ditentukan berdasarkan paparan alami
Bias	Lebih rendah karena randomisasi	Bisa lebih tinggi (cofounding mungkin terjadi)
Inferensi Kausal	Kuat	Relatif lebih lemah
Biaya dan Waktu	Mahal dan lama	Lebih murah dan lebih cepat
Contoh Studi	Uji klinis obat	Studi risiko merokok terhadap kanker

Dan beberapa perbandingan untuk penggunaan studi kohort dan studi kasus-kontrol akan dijelaskan dalam tabel di bawah ini :

Tabel 3. Penggunaan metode studi kohort dan studi kasus-kontrol
Faktor	Studi Kohort	Studi Kasus-Kontrol
Dimulai Dari	Paparan (misalnya perokok vs bukan perokok)	Outcome (misalnya kanker vs sehat)
Arah Studi	Prospektif / retrospektif	Retrospektif
Ukuran Efek	Relative Risk	Odds Ratio
Cocok Untuk	Penyakit umum	Penyakit langka
Biaya dan Waktu	Mahal dan lama	Lebih cepat dan murah
Bias	Lebih rendah	Lebih tinggi
Kemampuan Inferensi Kausal	Lebih kuat	Lebih lemah

5 TABEL KONTINGENSI 2X2

Tabel kontingensi adalah tabel yang digunakan untuk menyajikan distribusi frekuensi dari dua atau lebih variabel kategori. Tabel kontingensi 2x2 adalah bentuk paling sederhana dari tabel kontingensi yang digunakan untuk menganalisis hubungan antara dua variabel kategori.

Struktur dari tabel kontingensi adalah sebagai berikut :

Tabel 4. Struktur Tabel Kontingensi 2x2
	Y1	Y2	Total
X1	n11	n12	n1.
X2	n21	n22	n2.
Total	n.1	n.2	n..

Dimana :

n11 = Nilai pengamatan untuk kategori X ke-1 dan kategori Y ke-1
n12 = Nilai pengamatan untuk kategori X ke-1 dan kategori Y ke-2
n21 = Nilai pengamatan untuk kategori X ke-2 dan kategori Y ke-1
n22 = Nilai pengamatan untuk kategori X ke-2 dan kategori Y ke-2
n1. = Total pengamatan untuk kategori X ke-1
n2. = Total pengamatan untuk kategori X ke-2
n.1 = Total pengamatan untuk kategori Y ke-1
n.2 = Total pengamatan untuk kategori Y ke-2
n.. = Total seluruh nilai pengamatan

Misalkan variabel X terdiri dari I kategori dan variabel Y terdiri dari J kategori. Maka, tabel kontingensi antara X dan Y adalah tabel kontingensi I x J.

5.1 DISTRIBUSI FREKUENSI DAN PROPORSI

5.1.1 FREKUENSI ABSOLUT

Jumlah pengamatan dalam setiap sel tabel.

5.1.2 FREKUENSI RELATIF

Menunjukkan proporsi relatif suatu sel terhadap total populasi.

5.1.3 DISTRIBUSI MARGINAL

Distribusi total di setiap baris atau kolom / distribusi total dari masing-masing variabel kategorik.

5.1.4 DISTRIBUSI KONDISIONAL

Probabilitas suatu kejadian berdasarkan syarat tertentu.

5.2 DISTRIBUSI PELUANG DALAM TABEL KONTINGENSI 2X2

Untuk setiap tabel kontingensi dua arah, bila variabel X dan Y merupakan variabel acak dan pemilihan sampel dilakukan secara acak, maka tabel kontingensi tersebut memiliki distribusi peluang.

Tabel 5. Tabel Distribusi Peluang Tabel Kontingensi 2x2
	Y1	Y2	Total
X1	𝜋11	𝜋12	𝜋1.
X2	𝜋21	𝜋22	𝜋2.
Total	𝜋.1	𝜋.2	𝜋..

Dimana :

𝜋11 = Peluang bersama untuk kategori X ke-1 dan kategori Y ke-1
𝜋12 = Peluang bersama untuk kategori X ke-1 dan kategori Y ke-2
𝜋21 = Peluang bersama untuk kategori X ke-2 dan kategori Y ke-1
𝜋22 = Peluang bersama untuk kategori X ke-2 dan kategori Y ke-2
𝜋1. = Total peluang bersama untuk kategori X ke-1
𝜋2. = Total peluang bersama untuk kategori X ke-2
𝜋.1 = Total peluang bersama untuk kategori Y ke-1
𝜋.2 = Total peluang bersama untuk kategori Y ke-2
𝜋.. = 1

5.2.1 PELUANG BERSAMA (JOINT PROBABILITY)

Peluang bersama adalah probabilitas bahwa kedua variabel, X dan Y, terjadi secara bersamaan dalam suatu sel tabel kontingensi.

Peluang bersama dapat dihitung melalui rumus di bawah ini :

\[ P(X\cap Y) = \frac{n(I \cap J)}{n(S)} \]

5.2.2 PELUANG MARGINAL (MARGINAL PROBABILITY)

Peluang marginal adalah probabilitas suatu kejadian tunggal, baik itu kejadian kategori X atau Y, tanpa mempertimbangkan kejadian dari variabel lainnya.

Peluang marginal dapat dihitung melalui rumus di bawah ini :

\[ P(X) = \frac{n(I)}{n(S)} \]

5.2.3 PELUANG BERSYARAT (CONDITIONAL PROBABILITY)

Peluang bersyarat adalah probabilitas suatu kejadian terjadi dengan syarat kejadian lain telah terjadi. Langkah pertama yang harus dilakukan adalah menentukan variabel yang menjadi syarat terlebih dahulu.

Peluang bersyarat dapat dihitung melalui rumus di bawah ini :

\[ P(J \mid I) = \frac{P(I \cap J)}{P(I)} \]

Contoh soal menggunakan bantuan software R studio :

Seorang dokter melakukan pengujian obat sakit kepala kepada pasien berjumlah 300 orang. Obat yang digunakan yaitu Bodrex dan Paramex dimana masing-masing obat diberikan kepada 150 pasien. Setelah beberapa jam kemudian, diperiksa hasil kinerja dari obat-obat tersebut yang dilampirkan dalam tabel kontingensi berikut :

library(kableExtra)

# Dataset contoh
tabel <- matrix(c(100, 50, 80, 70), nrow = 2, byrow = TRUE,
               dimnames = list(Grup = c("Bodrex", "Paramex"),
                               Pusing = c("Ya", "Tidak")))
Total_tabel <- rbind(tabel, colSums(tabel))
Total_tabel <- cbind(Total_tabel, rowSums(Total_tabel))
dimnames(Total_tabel) <- list(c("Bodrex", "Paramex", "Total"), c("Ya", "Tidak", "Total"))

# Membuat tabel
kable(Total_tabel, caption = "Tabel Kontingensi Pengujian Obat", 
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE) %>%
  row_spec(nrow(Total_tabel), bold = TRUE, italic = FALSE)%>%
  collapse_rows(columns=1, valign="middle")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Kontingensi Pengujian Obat
	Ya	Tidak	Total
Bodrex	100	50	150
Paramex	80	70	150
Total	180	120	300

Peluang Bersama

n <- sum(tabel)

# Peluang Bersama
P_joint <- tabel / n
P_joint_df <- as.data.frame(P_joint)
P_joint_df$Grup <- rownames(P_joint_df)
P_joint_df <- P_joint_df[, c("Grup", "Ya", "Tidak")]
rownames(P_joint_df) <- NULL
P_joint_df %>%
  kable(caption = "Tabel Peluang Bersama", booktabs = TRUE, digits = 4) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Peluang Bersama
Grup	Ya	Tidak
Bodrex	0.3333	0.1667
Paramex	0.2667	0.2333

Melalui distribusi peluang bersama, dapat disimpulkan beberapa hal yakni :

Peluang dari pasien yang mengonsumsi Bodrex dan mengalami pusing adalah sebesar 0,333. Sementara yang mengonsumsi Paramex dan mengalami pusing hanya sebesar 0,267.
Peluang dari pasien yang mengonsumsi Bodrex dan tidak mengalami pusing hanya sebesar 0,167. Sementara yang mengonsumsi Paramex dan tidak mengalami pusing adalah sebesar 0,234.
Ini menunjukkan bahwa Paramex memiliki efektivitas penyembuhan yang lebih tinggi dibandingkan Bodrex berdasarkan hasil uji.

Peluang Marginal

P_marginal_grup <- rowSums(tabel) / n
P_marginal_pusing <- colSums(tabel) / n
P_marginal_df <- data.frame(
  Kategori = c("Bodrex", "Paramex", "Ya", "Tidak"),
  Peluang_Marginal = c(P_marginal_grup, P_marginal_pusing)
)

# Peluang Marginal
rownames(P_marginal_df) <- NULL
P_marginal_df %>%
  kable(caption = "Tabel Peluang Marginal", booktabs = TRUE, digits = 4) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Peluang Marginal
Kategori	Peluang_Marginal
Bodrex	0.5
Paramex	0.5
Ya	0.6
Tidak	0.4

Melalui distribusi peluang marginal, dapat disimpulkan beberapa hal yakni :

Peluang dari pasien yang mengonsumsi Bodrex dan Paramex adalah seimbang yaitu sama-sama bernilai 0,5.
Peluang marginal dari pasien yang mengalami gejala pusing setelah mengonsumi obat adalah sebesar 0,6, sedangkan yang tidak mengalami gejala pusing adalah sebesar 0,4. Ini menunjukkan bahwa secara keseluruhan, obat-obatan yang diuji cenderung tidak efektif bagi mayoritas pasien.

Peluang Bersyarat

# Peluang bersyarat (Pusing | Grup)
P_cond_baris <- tabel / rowSums(tabel)
P_cond_PG <- as.data.frame(P_cond_baris)
P_cond_PG$Grup <- rownames(P_cond_PG)
P_cond_PG <- P_cond_PG[, c("Grup", "Ya", "Tidak")]

# Peluang Bersyarat 1
rownames(P_cond_PG) <- NULL
P_cond_PG %>%
  kable(caption = "Tabel Peluang Bersyarat (Pusing | Grup)", booktabs = TRUE, digits = 4) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Peluang Bersyarat (Pusing | Grup)
Grup	Ya	Tidak
Bodrex	0.6667	0.3333
Paramex	0.5333	0.4667

# Peluang bersyarat (Grup | Pusing)
P_cond_kolom <- prop.table(tabel, margin = 2)
P_cond_GP <- as.data.frame(P_cond_kolom)
P_cond_GP$Grup <- rownames(P_cond_GP)
P_cond_GP <- P_cond_GP[, c("Grup", "Ya", "Tidak")]

# Peluang Bersyarat 2
rownames(P_cond_GP) <- NULL
P_cond_GP %>%
  kable(caption = "Tabel Peluang Bersyarat (Grup | Pusing)",
        booktabs = TRUE, digits = 4) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Peluang Bersyarat (Grup | Pusing)
Grup	Ya	Tidak
Bodrex	0.5556	0.4167
Paramex	0.4444	0.5833

Melalui distribusi peluang bersyarat, dapat disimpulkan beberapa hal yakni :

Peluang dari pasien yang mengalami gejala pusing menurut jenis obat yang dikonsumsi adalah sebesar 0,667 untuk pengguna Bodrex dan sebesar 0,533 untuk pengguna Paramex.
Terlihat dari pasien pengguna Bodrex menurut adanya gejala pusing adalah dominan mengalami. Sementara peluang dari pasien pengguna Paramex menurut adanya gejala pusing adalah dominan tidak mengalami.
Ini menunjukkan bahwa Paramex memiliki efektivitas penyembuhan yang lebih tinggi dibandingkan Bodrex berdasarkan hasil uji.

5.3 UKURAN ASOSIASI DALAM TABEL KONTINGENSI 2X2

Ukuran asosiasi adalah indikator statistik yang digunakan untuk mengukur kekuatan dan arah hubungan antara dua variabel, khususnya dalam konteks variabel kategorik seperti dalam tabel kontingensi yang memiliki lebih dari sama dengan 2 kategori.

Analisis asosiasi ini sering diterapkan dalam berbagai bidang seperti epidemiologi, riset sosial, dan eksperimen klinis seperti di bawah ini :

Epidemiologi : Meneliti hubungan antara kebiasaan merokok dan kanker paru-paru.
Eksperimen klinis : Mengukur efektivitas suatu pengobatan terhadap penyakit.
Riset sosial : Memeriksa hubungan antara tingkat pendidikan dan tingkat pekerjaan.

Beberapa jenis ukuran asosiasi dan fungsinya akan dijelaskan dalam tabel di bawah ini :

Tabel 6. Tabel Jenis Ukuran Asosiasi
Ukuran	Fungsi	Penelitian
Risk Difference	Mengukur selisih absolut risiko antar kelompok	Studi kohort & Eksperimen
Relative Risk	Perbandingan risiko relatif antara kelompok	Studi kohort & Eksperimen
Odds Ratio	Perbandingan odds kejadian antar kelompok	Studi kasus-kontrol

5.3.1 RISK DIFFERENCE (PERBEDAAN RISIKO)

Risk Difference menunjukkan selisih risiko / probabilitas antara kelompok yang terpapar dan tidak terpapar terhadap suatu kejadian.

Rumus Risk Difference dapat dituliskan dalam formula berikut :

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

a = individu terpapar dan mengalami outcome

b = individu terpapar dan tidak mengalami outcome

c = individu tidak terpapar dan mengalami outcome

d = individu tidak terpapar dan tidak mengalami outcome

Dan kriteria uji dalam Risk Difference adalah :

Jika RD > 0, maka risiko kejadian lebih tinggi di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika RD < 0, maka risiko kejadian lebih rendah di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika RD = 0, maka tidak ada asosiasi antara dua kelompok.

5.3.2 RELATIVE RISK (RISIKO RELATIF)

Relative Risk mengukur seberapa besar risiko kejadian pada kelompok terpapar relatif terhadap kelompok tidak terpapar.

Rumus Relative Risk dapat dituliskan dalam formula berikut :

\[ RR =\frac {\frac{a}{a+b}} {\frac{c}{c+d}} \]Dan kriteria uji dalam Relative Risk adalah :

Jika RR > 1, maka risiko kejadian lebih tinggi di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika RR < 1, maka risiko kejadian lebih rendah di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika RR = 1, maka tidak ada asosiasi antara dua kelompok.

5.3.3 ODDS RATIO (RASIO ODDS)

Odds Ratio menunjukkan perbandingan peluang kejadian antar kelompok.

Rumus Odds Ratio dapat dituliskan dalam formula berikut :

\[ OR = \frac{a \times d}{b \times c} \] Dan kriteria uji dalam Odds Ratio adalah :

Jika OR > 1, maka peluang kejadian lebih tinggi di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika OR < 1, maka peluang kejadian lebih rendah di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika OR = 1, maka tidak ada asosiasi antara dua kelompok.

Contoh soal menggunakan bantuan software R studio :

Sebuah survei dilakukan terhadap 200 orang untuk mengetahui hubungan antara pola hidup sehat dan hasil tes kesehatan (positif atau negatif). Berikut adalah data yang diperoleh :

# Dataset contoh
data <- matrix(c(65, 40, 20, 75), nrow=2, byrow=TRUE, dimnames=list(Pola_hidup=c("Buruk", "Baik"), Hasil_Tes=c("Negatif", "Positif")))  

# Buat tabel kontingensi
table_contingency <- as.table(data)
table_contingency

##           Hasil_Tes
## Pola_hidup Negatif Positif
##      Buruk      65      40
##      Baik       20      75

Risk Difference

# RD
RD <- function(table_contingency) {
  a <- table_contingency[1,1]
  b <- table_contingency[1,2]
  c <- table_contingency[2,1]
  d <- table_contingency[2,2]
  (a/(a+b))-(c/(c+d))
}
RD(table_contingency)

## [1] 0.4085213

Berdasarkan hasil perhitungan Risk Difference tersebut, diperoleh nilai RD sebesar 0,408. Karena nilai RD > 0, maka kelompok dengan pola hidup buruk memiliki risiko hasil tes negatif yang lebih tinggi dibandingkan kelompok dengan pola hidup baik .

Relative Risk

# RR
RR <- function(table_contingency) {
  a <- table_contingency[1,1]
  b <- table_contingency[1,2]
  c <- table_contingency[2,1]
  d <- table_contingency[2,2]
  (a/(a+b))/(c/(c+d))
}
RR(table_contingency)

## [1] 2.940476

Berdasarkan hasil perhitungan Relative Risk tersebut, diperoleh nilai RR sebesar 2,94. Karena nilai RR > 0, maka risiko hasil kesehatan negatif 2,94 kali lebih tinggi terjadi pada kelompok dengan pola hidup buruk dibandingkan kelompok dengan pola hidup baik .

Odds Ratio

# OR
OR <- function(table_contingency) {
  a <- table_contingency[1,1]
  b <- table_contingency[1,2]
  c <- table_contingency[2,1]
  d <- table_contingency[2,2]
  (a*d)/(b*c)
}
OR(table_contingency)

## [1] 6.09375

Berdasarkan hasil perhitungan Odds Ratio tersebut, diperoleh nilai OR sebesar 6,093. Karena nilai OR > 0, maka kelompok dengan pola hidup buruk memiliki peluang lebih tinggi untuk mendapatkan hasil tes negatif dibandingkan kelompok dengan pola hidup baik .

6 INFERENSI TABEL KONTINGENSI 2X2

Inferensi statistik adalah proses menarik kesimpulan tentang populasi berdasarkan data sampel. Dalam tabel kontingensi, inferensi digunakan untuk :

Mengestimasi ukuran asosiasi antara 2 variabel.
Mengujinya untuk melihat apakah asosiasi itu signifikan secara statistik.

6.1 ESTIMASI

Estimasi digunakan untuk mengukur besarnya asosiasi antara 2 variabel. Estimasi terbagi menjadi :

6.1.1 ESTIMASI TITIK

Estimasi titik menunjukkan nilai tunggal terbaik dari parameter populasi berdasarkan data sampel.

Rumus estimasi titik dapat dituliskan dalam formula berikut :

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

p topi = estimasi titik proporsi

x = jumlah individu dalam kategori tertentu

n = total jumah individu dalam sampel

6.1.2 ESTIMASI INTERVAL

Estimasi juga biasanya disertai confidence interval (CI) untuk menggambarkan ketidakpastian nilai dalam estimasi yang dipercaya mengandung parameter populasi dengan tingkat kepercayaan tertentu.

Rumus estimasi interval dapat dituliskan dalam formula berikut :

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

z = nilai distribusi normal standar untuk tingkat kepercayaan tertentu

p topi = estimasi titik proporsi

n = total jumah individu dalam sampel

6.2 UJI HIPOTESIS

Uji hipotesis adalah prosedur statistik untuk mengevaluasi klaim atau pernyataan tentang parameter populasi berdasarkan data sampel. Tujuan utamanya adalah menentukan apakah terdapat cukup bukti dari sampel untuk menolak H0 pada tingkat kepercayaan tertentu.

6.2.1 UJI PROPORSI

Uji proporsi digunakan untuk menguji apakah proporsi suatu kejadian dalam populasi sama dengan nilai tertentu atau berbeda antar dua kelompok. Dalam tabel kontingensi 2x2, akan digunakan uji Z dua proporsi.

Dan hipotesis untuk uji proporsi adalah :

H0 = Tidak terdapat perbedaan proporsi antara dua kelompok.
H1 = Terdapat perbedaan proporsi antara dua kelompok.

Statistik uji proporsi dua sampel dapat dituliskan dalam formula berikut : \[ Z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{\hat{p}(1 - \hat{p}) \left( \frac{1}{n_1} + \frac{1}{n_2} \right)}} \]

Dan kriteria uji dalam uji proporsi adalah :

Jika p-value < 0,05 atau z-value > 0,05 maka tolak H0 dan dapat disimpulkan bahwa terdapat perbedaan proporsi antara dua kelompok.
Jika p-value > 0,05 atau z-value < 0,05 maka terima H0 dan dapat disimpulkan bahwa tidak terdapat perbedaan proporsi antara dua kelompok.

6.2.2 UJI ASOSIASI

Uji asosiasi digunakan untuk melihat apakah terdapat asosiasi antara dua variabel kategori.

Tiga ukuran umum dari tabel kontingensi 2x2 di antaranya :

Risk Difference (RD)
Relative Risk (RR)
Odds Ratio (OR)

Dan hipotesis untuk uji asosiasi adalah :

H0 = Tidak terdapat asosiasi antara dua kelompok.
H1 = Terdapat perbedaan asosiasi antara dua kelompok.

Statistik uji Z untuk Risk Difference dapat dituliskan dalam formula berikut :

\[ Z_{RD}= {\frac{\frac{a}{a+b}-\frac{c}{c+d}}{\sqrt{ \frac{ \hat{p}_1 (1 - \hat{p}_1) }{n_{1\cdot}} + \frac{ \hat{p}_2 (1 - \hat{p}_2) }{n_{2\cdot}} }}} \]

Dan kriteria uji dalam uji asosiasi adalah :

Jika p-value < 0,05 atau z-value > 0,05 maka tolak H0 dan dapat disimpulkan bahwa terdapat asosiasi antara dua kelompok.
Jika p-value > 0,05 atau z-value < 0,05 maka terima H0 dan dapat disimpulkan bahwa tidak terdapat asosiasi antara dua kelompok.

Contoh soal menggunakan bantuan software R studio :

Sebuah sekolah ingin mengetahui apakah ada perbedaan proporsi kelulusan ujian matematika nasional antara siswa kelas reguler dan kelas akselerasi. Diperoleh data sebagai berikut :

# Dataset contoh
tabel <- matrix(c(42, 18, 38, 2), nrow = 2, byrow = TRUE,
               dimnames = list(Kelas = c("Reguler", "Akselerasi"),
                               Lulus = c("Ya", "Tidak")))
Total_tabel <- rbind(tabel, colSums(tabel))
Total_tabel <- cbind(Total_tabel, rowSums(Total_tabel))
dimnames(Total_tabel) <- list(c("Reguler", "Akselerasi", "Total"), c("Lulus", "Tidak Lulus", "Total"))

# Membuat tabel
kable(Total_tabel, caption = "Tabel Kontingensi Kelulusan Siswa", 
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE) %>%
  row_spec(nrow(Total_tabel), bold = TRUE, italic = FALSE)%>%
  collapse_rows(columns=1, valign="middle")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Kontingensi Kelulusan Siswa
	Lulus	Tidak Lulus	Total
Reguler	42	18	60
Akselerasi	38	2	40
Total	80	20	100

# Uji Proporsi
prop_test <- prop.test(x = c(42, 38), n = c(60, 40), correct = FALSE)
print(prop_test)

## 
##  2-sample test for equality of proportions without continuity correction
## 
## data:  c(42, 38) out of c(60, 40)
## X-squared = 9.375, df = 1, p-value = 0.0022
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.3841896 -0.1158104
## sample estimates:
## prop 1 prop 2 
##   0.70   0.95

Berdasarkan hasil uji proporsi di atas, diperoleh p-value senilai 0,0022 ( < 0,05) sehingga H0 ditolak dan dapat disimpulkan bahwa terdapat perbedaan proporsi antara siswa kelas reguler dan kelas akselerasi.

# Uji Asosiasi
# Contoh dataset
a <- 42
b <- 18
c <- 38
d <- 2

# Risk Difference
p1 <- (a/(a+b))
p2 <- (c/(c+d))
rd <- p1 - p2
se_rd <- sqrt((p1*(1-p1)/(a+b))+p2*((1-p2)/(c+d)))
z_rd <- rd/se_rd

# Relative Risk
rr <- (a/(a+b)) / (c/(c+d))
se_ln_rr <- sqrt((1/a) - (1/(a+b)) + (1/c) - (1/(c+d)))
z_rr <- log(rr) / se_ln_rr

# Odds Ratio
or <- (a * d) / (b * c)
se_ln_or <- sqrt((1/a) + (1/b) + (1/c) + (1/d))
z_or <- log(or) / se_ln_or

#Hasil
hasil <- list(
  Z_RD = z_rd,
  Z_RR = z_rr,
  Z_OR = z_or
)

hasil_df <- data.frame(
  Ukuran = names(hasil),
  Nilai = unlist(hasil)
)

library(kableExtra)
rownames(hasil_df) <- NULL
hasil_df %>%
  kable(
    caption = "Tabel Ukuran Asosiasi",
    digits = 4,
    align = "lc",
    booktabs = TRUE
  ) %>%
  kable_styling(full_width = FALSE, position = "center")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Ukuran Asosiasi
Ukuran	Nilai
Z_RD	-3.6515
Z_RR	-3.3204
Z_OR	-2.6947

chi_test <- chisq.test(tabel)
print(chi_test)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  tabel
## X-squared = 7.8776, df = 1, p-value = 0.005005

Berdasarkan hasil uji asosiasi di atas, diperoleh p-value senilai 0,005 ( < 0,05) sehingga H0 ditolak dan dapat disimpulkan bahwa jenis kelas mempengaruhi kelulusan secara signifikan.

Berdasarkan hasil perhitungan statistik uji Risk Difference tersebut, diperoleh nilai RD sebesar -3,65. Karena nilai RD < 0, maka kelompok siswa yang berada di kelas akselerasi memiliki risiko lulus tes ujian matematika yang lebih tinggi dibandingkan kelompok siswa yang berada di kelas reguler .

Berdasarkan hasil perhitungan Relative Risk tersebut, diperoleh nilai RR sebesar -3,32. Karena nilai RR < 0, maka kelulusan ujian matematika 3,32 kali lebih tinggi teradi pada kelompok siswa yang berada di kelas akselerasi dibandingkan kelompok siswa yang berada di kelas reguler.

Berdasarkan hasil perhitungan Odds Ratio tersebut, diperoleh nilai OR sebesar -2,69. Karena nilai OR < 0, maka kelompok siswa yang berada di kelas akselerasi memiliki peluang lebih tinggi untuk lulus ujian matematika dibandingkan kelompok siswa yang berada di kelas reguler.

6.2.3 UJI INDEPENDENSI

Uji independensi digunakan untuk menetukan apakah terdapat asosiasi antara dua variabel kategori. Dalam tabel kontingens, kita menguji apakah distribusi suatu variabel tergantung pada kategori dari variabel lainnya.

6.2.3.1 UJI CHI-SQUARE

Uji chi-square digunakan untuk menguji apakah ada hubungan antara dua variabel kategori.

Hipotesis dalam uji chi-square adalah :

H0 = Tidak ada hubungan antara variabel kategori
H1 = Terdapat hubungan antara variabel kategori

Dan kriteria uji dalam uji chi-square adalah :

Jika p-value < 0,05 maka tolak H0 dan dapat disimpulkan bahwa terdapat hubungan antara variabel kategori.
Jika p-value > 0,05 maka terima H0 dan dapat disimpulkan bahwa tidak terdapat hubungan yang signifikan antara variabel kategori.

Rumus Chi-Square dapat dituliskan dalam formula berikut : \[ \chi^2 = \sum \frac{(O - E)^2}{E} \]

dimana :

O = nilai pengamatan
E = nilai yang diharapkan dan dapat diperoleh melalui formula :

\[ E = \frac{(Jumlah Baris × Jumlah Kolom )}{Jumlah Keseluruhan} \]

Partisi Chi-Square

Partisi chi-square adalah teknik untuk memecah (mempartisi) nilai chi-square total dari suatu tabel kontingensi besar menjadi beberapa komponen chi-square yang lebih kecil. Tujuan utamanya adalah untuk mengidentifikasi sumber spesifik yang menjadi pengaruh signifikan antara kategori pengamatan di dalam tabel kontingensi.
Uji Likelihood Ratio

Uji Likelihood Ratio merupakan alternatif dari uji chi-square dan lebih baik digunakan jika data kecil atau tidak proporsional.

Rumus Likelihood Ratio dapat dituliskan dalam formula berikut : \[ G^2 =2 \sum O_{ij}log(\frac{O_{ij}}{E_{ij}}) \]

dimana :

Oij = frekuensi observasi dalam tabel kontingensi

Eij = frekuensi ekspektasi yang diperoleh melalui rumus :

\[ E_{ij}=n(p_{i+}+p_{+j}) \]

Contoh soal menggunakan bantuan software R studio :

Sebuah survei dilakukan terhadap 200 orang untuk mengetahui hubungan antara usia dan hasil tes kesehatan diare (positif atau negatif). Berikut adalah data yang diperoleh :

# Dataset contoh
data <- matrix(c(54,5,17,66,24,34), nrow=2, byrow=TRUE)
colnames(data) <- c("Anak-anak", "Remaja", "Dewasa")
rownames (data) <- c("Positif", "Negatif")

# Buat tabel kontingensi
table_contingency <- as.table(data)
print(addmargins(table_contingency))

##         Anak-anak Remaja Dewasa Sum
## Positif        54      5     17  76
## Negatif        66     24     34 124
## Sum           120     29     51 200

Uji Partisi Chi-Square

# Uji Chi-Square
chi_test <- chisq.test(table_contingency)
print(chi_test)

## 
##  Pearson's Chi-squared test
## 
## data:  table_contingency
## X-squared = 8.2714, df = 2, p-value = 0.01599

# Uji Partisi
# Data observasi 1
data_matrix <- matrix(c(54,5,66,24), nrow=2, byrow=TRUE)
colnames(data_matrix) <- c("Anak-anak", "Remaja")
rownames (data_matrix) <- c("Positif", "Negatif")
data_matrix

##         Anak-anak Remaja
## Positif        54      5
## Negatif        66     24

chi_test1 <- chisq.test(data_matrix)
print (chi_test1)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_matrix
## X-squared = 6.4085, df = 1, p-value = 0.01136

# Data observasi 2
data_matrix2 <- matrix(c(54,17,66,34), nrow=2, byrow=TRUE)
colnames(data_matrix2) <- c("Anak-anak", "Dewasa")
rownames (data_matrix2) <- c("Positif", "Negatif")
data_matrix2

##         Anak-anak Dewasa
## Positif        54     17
## Negatif        66     34

chi_test2 <- chisq.test(data_matrix2)
print (chi_test2)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_matrix2
## X-squared = 1.5545, df = 1, p-value = 0.2125

# Data observasi 3
data_matrix3 <- matrix(c(5,17,24,34), nrow=2, byrow=TRUE)
colnames(data_matrix3) <- c("Remaja", "Dewasa")
rownames (data_matrix3) <- c("Positif", "Negatif")
data_matrix3

##         Remaja Dewasa
## Positif      5     17
## Negatif     24     34

chi_test3 <- chisq.test(data_matrix3)
print (chi_test3)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_matrix3
## X-squared = 1.6619, df = 1, p-value = 0.1973

# Data observasi 4
data_matrix4 <- matrix(c(54,22,66,58), nrow=2, byrow=TRUE)
colnames(data_matrix4) <- c("Anak-anak", "Remaja+Dewasa")
rownames (data_matrix4) <- c("Positif", "Negatif")
data_matrix4

##         Anak-anak Remaja+Dewasa
## Positif        54            22
## Negatif        66            58

chi_test4 <- chisq.test(data_matrix4)
print (chi_test4)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_matrix4
## X-squared = 5.5187, df = 1, p-value = 0.01881

#Hasil
list(Chi_Square_Partisi1 = chi_test1, Chi_Square_Partisi2 = chi_test2, Chi_Square_Partisi3 = chi_test3,Chi_Square_Partisi4 = chi_test4)

## $Chi_Square_Partisi1
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_matrix
## X-squared = 6.4085, df = 1, p-value = 0.01136
## 
## 
## $Chi_Square_Partisi2
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_matrix2
## X-squared = 1.5545, df = 1, p-value = 0.2125
## 
## 
## $Chi_Square_Partisi3
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_matrix3
## X-squared = 1.6619, df = 1, p-value = 0.1973
## 
## 
## $Chi_Square_Partisi4
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_matrix4
## X-squared = 5.5187, df = 1, p-value = 0.01881

Berdasarkan hasil perhitungan uji chi-square tersebut, diperoleh nilai p-value senilai 0,015. Karena nilai p-value < 0,05 maka dapat disimpulkan bahwa H0 ditolak atau terdapat hubungan antara usia dan hasil tes kesehatan.

Sementara untuk uji partisi chi-square, diperoleh beberapa hasil sebagai berikut :

Uji anak-anak vs remaja menunjukkan hasil signifikan. ( p-value = 0,113 < alpha = 0,05)
Uji anak-anak vs dewasa menunjukkan hasil tidak signifikan. ( p-value = 0,212 > alpha = 0,05)
Uji remaja vs dewasa menunjukkan hasil tidak signifikan. ( p-value = 0,197 > alpha = 0,05)
Uji anak-anak vs remaja+dewasa menunjukkan hasil signifikan. ( p-value = 0,018 < alpha = 0,05)

Dengan mempartisi chi-square secara satu per satu, diketahui bahwa sumber spesifik yang menjadi pengaruh signifikan diantara kategori usia adalah anak-anak.

Likelihood Ratio

# Likelihood Ratio
data_expected <- chisq.test(data)$expected

G2 <- 2*sum(data*log(data/data_expected))
crit_value <- qchisq(0.95, df=1)

#Hasil
list(G2 = G2, Critical_Value = crit_value, Decision = ifelse(G2 > crit_value, "Reject H0", "Failed to Reject H0"))

## $G2
## [1] 8.885695
## 
## $Critical_Value
## [1] 3.841459
## 
## $Decision
## [1] "Reject H0"

Berdasarkan hasil perhitungan uji likelihood ratio tersebut, diperoleh nilai hitung sebesar 8,88. Karena nilai hitung > 3,84 maka dapat disimpulkan bahwa H0 ditolak atau terdapat hubungan antara usia dan hasil tes kesehatan.

6.2.3.2 UJI FISHER’S EXACT

Fisher’s Exact Test adalah uji statistik yang digunakan untuk menguji asosiasi antara dua variabel kategori, terutama ketika jumlah sampel kecil dan asumsi uji Chi-square tidak terpenuhi.

Hipotesis dalam uji Fisher’s Exact adalah :

H0 = Tidak ada hubungan antara variabel kategori
H1 = Terdapat hubungan antara variabel kategori

Dan kriteria uji dalam uji Fisher’s Exact adalah :

Jika p-value < 0,05 maka tolak H0 dan dapat disimpulkan bahwa terdapat hubungan antara variabel kategori.
Jika p-value > 0,05 maka terima H0 dan dapat disimpulkan bahwa tidak terdapat hubungan yang signifikan antara variabel kategori.

Rumus Fisher’s Exact dapat dituliskan dalam formula berikut : \[ P = \frac{(a+b)!(c+d)!(a+c)!(b+d)!}{a!b!c!d!n!} \]

Keunggulan
- Tidak bergantung pada asumsi distribusi Chi-square atau normalitas.
- Sangat akurat untuk data dengan frekuensi kecil.
- Cocok untuk sampel kecil bahkan ketika total hanya 4 atau 5 observasi.
Kekurangan
- Komputasi berat untuk tabel kontingensi yang besar.
- Bisa menghasilkan p-value yang lebih besar daripada seharusnya.

Contoh soal menggunakan bantuan software R studio :

Seorang guru ingin mengetahui apakah metode belajar mempengaruhi hasil ujian siswa. Sebanyak 20 siswa dibagi menjadi dua kelompok dengan perolehan data sebagai berikut :

# Data
data <- matrix(c(8, 2, 4, 6), nrow = 2, byrow = TRUE)
colnames(data) <- c("Lulus", "Tidak Lulus")
rownames(data) <- c("Interaktif", "Ceramah")

# Uji Fisher
fisher.test(data)

## 
##  Fisher's Exact Test for Count Data
## 
## data:  data
## p-value = 0.1698
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
##   0.6026805 79.8309210
## sample estimates:
## odds ratio 
##   5.430473

Berdasarkan hasil perhitungan uji Fisher’s Exact tersebut, diperoleh nilai p-value sebesar 0,169. Karena p-value > 0,05 maka dapat disimpulkan bahwa H0 diterima atau tidak terdapat hubungan antara metode belajar dan hasil ujian siswa.

6.3 ANALISIS RESIDU

Analisis residu dilakukan setelah uji chi-square untuk mengetahui sel mana dalam tabel kontingensi yang berkontribusi besar terhadap hubungan antara variabel kategori. Jika suatu sel memiliki residual yang kecil (positif atau negatif), maka nilai observasi tersebut mendekati nilai ekspektasi sehingga sel tersebut tidak banyak menyumbang untuk hubungan anatra variabel.

Berikut adalah beberapa catatan mengenai kesimpulan asosiasi menurut analisis residualnya :

Jika residual kecil , tidak ada perbedaan signifikan antara jumlah pengamatan dengan jumlah yang diprediksi sehingga antar variabel tidak memiliki hubungan kuat.
Jika residual positif besar, jumlah pengamatan jauh lebih tinggi dari yang diprediksi sehingga antar variabel memiliki hubungan positif yang kuat.
Jika residual negatif besar, jumlah pengamatan jauh lebih rendah dari yang diprediksi sehingga antar variabel memiliki hubungan negatif atau tidak ada asosiasi.
Jika residual kecil (0), maka tidak ada asosiasi antar variabel.

6.3.1 JENIS RESIDUAL DAN DETEKSI OUTLIER

6.3.1.1 PEARSON RESIDUAL (RESIDU STANDAR)

Pearson residual dapat dihitung melalui rumus berikut :

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

dimana :

Oij = nilai pengamatan pada baris ke-i kolom ke-j

Eij = nilai harapan pada baris ke-i kolom ke-j

6.3.1.2 STANDARDIZED RESIDUAL (Z-SCORE RESIDUAL)

Standardized residual dapat dihitung melalui rumus berikut :

\[ r_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}(1 - p_{i+})(1 - p_{+j})}} \]

dimana :

pi+ / p+j = probabilitas marginal dari baris dan kolom

6.3.1.3 DETEKSI OUTLIER

Deteksi outlier bertujuan untuk mencari sel-sel yang memiliki nilai penyimpangan paling ekstrem antara frekuensi observasi dan frekuensi harapan. Dalam analisis residual, outlier bukan berarti angka tunggal yang jauh dari rata-rata, melainkan sel pada tabel kontingensi yang berkontribusi besar terhadap nilai chi-square karena residunya yang sangat besar.

Dan kriteria uji dalam analisis residual adalah :

Standardized Residual (SR) ≥ |2|, maka terdapat perbedaan yang signifikan.
Standardized Residual (SR) ≥ |3|, maka data merupakan outlier.

Contoh soal menggunakan bantuan software R studio :

Sebuah survei dilakukan terhadap 200 orang untuk mengetahui hubungan antara usia dan hasil tes kesehatan diare (positif atau negatif). Berikut adalah data yang diperoleh :

# Dataset contoh
data <- matrix(c(54,5,17,66,24,34), nrow=2, byrow=TRUE)
expected <- chisq.test(data)$expected

# Pearson Residual
pearson_residual <- (data - expected)/sqrt(expected)

# Standardized Residual
row_sum <- rowSums(data)
col_sum <- colSums (data)
total_sum <- sum(data)

# Probabilitas marginal
p_row <- row_sum / total_sum
p_col <- col_sum / total_sum

# Expected frequencies
expected <- outer(row_sum, col_sum) / total_sum

# Standardized residuals
standardized_residual <- matrix(NA, nrow = nrow(data), ncol = ncol(data))

for (i in 1:nrow(data)) {
  for (j in 1:ncol(data)) {
    pij <- expected[i, j]
    if (pij > 0) {
      denom <- sqrt(pij * (1 - p_row[i]) * (1 - p_col[j]))
      if (denom > 0) {
        standardized_residual[i, j] <- (data[i, j] - pij) / denom
      }
    }
  }
}

# Hasil
list(
Pearson_Residual = pearson_residual,
Standardized_Residual = standardized_residual
)

## $Pearson_Residual
##            [,1]      [,2]       [,3]
## [1,]  1.2439326 -1.813450 -0.5406299
## [2,] -0.9738517  1.419717  0.4232491
## 
## $Standardized_Residual
##           [,1]      [,2]       [,3]
## [1,]  2.497877 -2.490731 -0.7954742
## [2,] -2.497877  2.490731  0.7954742

# Deteksi Outlier
library(knitr)
library(kableExtra)

# Data Pearson Residual
pearson_residual <- matrix(c(
  1.2439326, -1.813450, -0.5406299,
  -0.9738517, 1.419717, 0.4232491
), nrow = 2, byrow = TRUE)
colnames(pearson_residual) <- c("Anak-anak", "Remaja", "Dewasa")
rownames(pearson_residual) <- c("Positif", "Negatif")

# Data Standardized Residual
standardized_residual <- matrix(c(
  2.497877, -2.490731, -0.7954742,
  -2.497877, 2.490731, 0.7954742
), nrow = 2, byrow = TRUE)
colnames(standardized_residual) <- c("Anak-anak", "Remaja", "Dewasa")
rownames(standardized_residual) <- c("Positif", "Negatif")

kable(pearson_residual, caption = "Tabel Pearson Residual") %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Pearson Residual
	Anak-anak	Remaja	Dewasa
Positif	1.2439326	-1.813450	-0.5406299
Negatif	-0.9738517	1.419717	0.4232491

kable(standardized_residual, caption = "Tabel Standardized Residual") %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Standardized Residual
	Anak-anak	Remaja	Dewasa
Positif	2.497877	-2.490731	-0.7954742
Negatif	-2.497877	2.490731	0.7954742

Berdasarkan hasil analisis residu tersebut, diperoleh beberapa hasil sebagai berikut :

Anak-anak & Positif memiliki residual negatif dan positif terbesar (-2,49 , 2,49). Karena nilai SR > ±1,96 maka terdapat hubungan yang kuat antara kategori usia dan hasil tes kesehatan. Ini menunjukkan bahwa kelompok anak-anak mengalami lebih sedikit kegagalan dari yang diprediksi oleh model.
Dan karena (SR) ≤ |3|, maka tidak terdapat data outlier.

7 TABEL KONTINGENSI TIGA ARAH

Tabel kontingensi tiga arah menyajikan frekuensi pengamatan dari tiga variabel kategori (X, Y, dan Z) yang bertujuan untuk menilai interaksi antara dua variabel dikondisikan pada level variabel ketiga secara simultan dan menguji independensi bersyarat serta pengaruh perancu (cofounder).

7.1 TABEL PARSIAL DAN MARGINAL

Tabel parsial merupakan tabel yang diperoleh dengan mengelompokkan dua variabel berdasarkan satu variabel yang dikendalikan, sedangkan tabel marginal adalah tabel yang mengabaikan variabel yang dikendalikan tersebut dengan cara menjumlahkannya.

Kegunaan tabel marginal adalah untuk melihat pola asosiasi secara agregat, tetapi tabel parsial biasanya lebih diutamakan untuk mendapatkan kesimpulan yang lebih akurat.

Contoh soal menggunakan bantuan software R studio :

Sebuah survei ingin mengkaji apakah preferensi minum kopi (Ya/Tidak) tergantung pada status pekerjaan (Mahasiswa/Pekerja) dengan mempertimbangkan jenis kelamin (Laki-laki/Perempuan). Berikut adalah data yang diperoleh :

# Dataset contoh
data <- array(c(30,20,25,25,35,15,20,30),
dim=c(2,2,2),
dimnames=list(
Status=c("Mahasiswa", "Pekerja"),
Gender=c("Laki-laki","Perempuan"),
Kopi=c("Ya","Tidak")
))
data

## , , Kopi = Ya
## 
##            Gender
## Status      Laki-laki Perempuan
##   Mahasiswa        30        25
##   Pekerja          20        25
## 
## , , Kopi = Tidak
## 
##            Gender
## Status      Laki-laki Perempuan
##   Mahasiswa        35        20
##   Pekerja          15        30

library(knitr)
library(kableExtra)

# Tabel parsial
freq_parsial_ya <- data[,,"Ya"]
freq_parsial_tidak <- data[,,"Tidak"]
kable(freq_parsial_ya,caption = "Tabel Parsial Suka Kopi",
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Parsial Suka Kopi
	Laki-laki	Perempuan
Mahasiswa	30	25
Pekerja	20	25

kable(freq_parsial_tidak,caption = "Tabel Parsial Tidak Suka Kopi",
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Parsial Tidak Suka Kopi
	Laki-laki	Perempuan
Mahasiswa	35	20
Pekerja	15	30

# Tabel marginal
freq_marginal_X <- apply(data,1,sum)
df_marginal_X <- data.frame(
  Status = names(freq_marginal_X),
  Jumlah = as.vector(freq_marginal_X)
)
freq_marginal_Z <- apply(data,3,sum)
names(freq_marginal_Z) <- c("Laki-laki","Perempuan")
df_marginal_Z <- data.frame(
  Gender = names(freq_marginal_Z),
  Jumlah = as.vector(freq_marginal_Z)
)

kable(df_marginal_X, caption = "Tabel Marginal berdasarkan Status",
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Marginal berdasarkan Status
Status	Jumlah
Mahasiswa	110
Pekerja	90

kable(df_marginal_Z, caption = "Tabel Marginal berdasarkan Gender",
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Marginal berdasarkan Gender
Gender	Jumlah
Laki-laki	100
Perempuan	100

7.2 DISTRIBUSI PELUANG

7.2.1 PELUANG BERSAMA (JOINT PROBABILITY)

Peluang bersama mengukur seberapa besar kemungkinan tiga kategori terjadi bersamaan dalam suatu sel tabel kontingensi.

Peluang bersama dapat dihitung melalui rumus di bawah ini :

\[ P(X = i, Y = j, Z = k) = \frac{n\_{ijk}}{N} \]

7.2.2 PELUANG MARGINAL (MARGINAL PROBABILITY)

Peluang marginal adalah peluang dari satu kejadian (variabel) tanpa mempertimbangkan variabel lain.

Peluang marginal dapat dihitung melalui rumus di bawah ini :

\[ P(Y = j \mid X = i, Z = k) = \frac{P(X = i, Y = j, Z = k)}{P(X = i, Z = k)} = \frac{n_{ijk}}{\sum_j n_{ijk}} \]

7.2.3 PELUANG BERSYARAT (CONDITIONAL PROBABILITY)

Peluang bersyarat adalah peluang suatu kejadian terjadi, dengan asumsi kejadian lain telah diketahui terjadi.

Peluang bersyarat dapat dihitung melalui rumus di bawah ini :

\[ P(Y = j \mid X = i, Z = k) = \frac{P(X = i, Y = j, Z = k)}{P(X = i, Z = k)} = \frac{n_{ijk}}{\sum_j n_{ijk}} \]

Contoh soal menggunakan bantuan software R studio :

Sebuah penelitian dilakukan untuk mengetahui hubungan antara jenis kelamin, usia, dan jenis penyakit (Diare/Gastritis). Penelitian ini melibatkan sejumlah pasien yang dikategorikan berdasarkan:

Jenis Kelamin: Laki-laki / Perempuan
Usia: Anak-anak / Remaja / Dewasa
Jenis Penyakit: Diare / Gastritis

Data survei yang diperoleh adalah sebagai berikut:

library(magrittr) 
library(kableExtra)
library(dplyr)
library(knitr)

# Dataset contoh
data <- data.frame(
  Jenis_Kelamin = rep(c("Laki-laki", "Perempuan"), each = 3),
  Usia = rep(c("Anak-anak", "Remaja", "Dewasa"), times = 2),
  Diare = c(50,45,35,20,25,30),
  Gastritis = c(5,10,15,40,35,50)
)

data$Total_Baris <- data$Diare + data$Gastritis
total_kolom <- colSums(data[, 3:5])
total_kolom_df <- data.frame(Jenis_Kelamin = "Total Kolom", Usia = "", t(total_kolom))

data <- rbind(data, total_kolom_df)
data[3:5] <- lapply(data[3:5], as.numeric)

kable(data, caption = "Tabel Kontingensi Penyakit berdasarkan Jenis Kelamin dan Usia", 
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE) %>%
  row_spec(nrow(data), bold = TRUE, italic = FALSE)%>%
  collapse_rows(columns=1, valign="middle")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Kontingensi Penyakit berdasarkan Jenis Kelamin dan Usia
Jenis_Kelamin	Usia	Diare	Gastritis	Total_Baris
Laki-laki	Anak-anak	50	5	55
	Remaja	45	10	55
	Dewasa	35	15	50
Perempuan	Anak-anak	20	40	60
	Remaja	25	35	60
	Dewasa	30	50	80
Total Kolom		205	155	360

Peluang Bersama

prob_bersama <- data %>%
  mutate(
    P_Diare = Diare / sum(Diare),
    P_Gastritis = Gastritis / sum(Gastritis)
  )
kable(prob_bersama, caption = "Peluang Bersama Tiap Kategori", booktabs = TRUE) %>%
kable_styling(full_width = FALSE) %>%
row_spec(nrow(data))%>%
collapse_rows(columns=1, valign="middle")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Peluang Bersama Tiap Kategori
Jenis_Kelamin	Usia	Diare	Gastritis	Total_Baris	P_Diare	P_Gastritis
Laki-laki	Anak-anak	50	5	55	0.1219512	0.0161290
	Remaja	45	10	55	0.1097561	0.0322581
	Dewasa	35	15	50	0.0853659	0.0483871
Perempuan	Anak-anak	20	40	60	0.0487805	0.1290323
	Remaja	25	35	60	0.0609756	0.1129032
	Dewasa	30	50	80	0.0731707	0.1612903
Total Kolom		205	155	360	0.5000000	0.5000000

Melalui distribusi peluang bersama, dapat disimpulkan beberapa hal yakni :

Peluang terbesar untuk penyakit diare kategori laki-laki adalah di masa anak-anak yaitu sebesar 0,12. Sedangkan, peluang terbesar untuk penyakit gastritis kategori laki-laki adalah di masa dewasa yaitu sebesar 0,048.
Peluang terbesar untuk penyakit diare kategori perempuan adalah di masa dewasa yaitu sebesar 0,073. Sedangkan, peluang terbesar untuk penyakit gastritis kategori perempuan adalah di masa dewasa yaitu sebesar 0,16.
Ini menunjukkan bahwa pola kejadian penyakit berbeda antara laki-laki dan perempuan serta bervariasi di setiap usia. Dengan kata lain, jenis kelamin dan usia merupakan faktor penting terhadap jenis penyakit.

Peluang Marginal

prob_marginal <- data %>%
  summarise(
    P_Diare = sum(Diare) / sum(Total_Baris),
    P_Gastritis = sum(Gastritis) / sum(Total_Baris)
  )

kable(prob_marginal, caption = "Peluang Marginal", booktabs = TRUE)%>%
  kable_styling(full_width = FALSE) %>%
  collapse_rows(columns=1, valign="middle")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Peluang Marginal
P_Diare	P_Gastritis
0.5694444	0.4305556

Melalui distribusi peluang marginal, dapat disimpulkan bahwa secara keseluruhan, jenis penyakit diare lebih banyak terjadi dibandingkan penyakit gastritis.

Peluang Bersyarat

total_per_gender <- data %>%
  group_by(Jenis_Kelamin) %>%
  summarise(
    total_diare = sum(Diare),
    total_gastritis = sum(Gastritis)
  )

prob_bersyarat <- data %>%
  left_join(total_per_gender, by = "Jenis_Kelamin") %>%
  mutate(
    P_Diare_Gender = round(Diare / total_diare, 3),
    P_Gastritis_Gender = round(Gastritis / total_gastritis, 3)
  )

kable(prob_bersyarat, caption = "Peluang Bersyarat per Gender", booktabs = TRUE) %>%
  kable_styling(full_width = FALSE) %>%
  collapse_rows(columns = 1, valign = "middle")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Peluang Bersyarat per Gender
Jenis_Kelamin	Usia	Diare	Gastritis	Total_Baris	total_diare	total_gastritis	P_Diare_Gender	P_Gastritis_Gender
Laki-laki	Anak-anak	50	5	55	130	30	0.385	0.167
	Remaja	45	10	55	130	30	0.346	0.333
	Dewasa	35	15	50	130	30	0.269	0.500
Perempuan	Anak-anak	20	40	60	75	125	0.267	0.320
	Remaja	25	35	60	75	125	0.333	0.280
	Dewasa	30	50	80	75	125	0.400	0.400
Total Kolom		205	155	360	205	155	1.000	1.000

Melalui distribusi peluang bersyarat, dapat disimpulkan beberapa hal yakni :

Untuk penyakit diare kategori laki-laki, peluang tertinggi adalah di masa anak-anak (0,385) dan menurun di masa dewasa (0,269).
Untuk penyakit gastritis kategori laki-laki, peluang tertinggi adalah di masa dewasa (0,5).
Untuk penyakit diare dan gastritis kategori perempuan, peluang tertinggi adalah di masa dewasa (0,4).

7.3 UKURAN ASOSIASI DALAM TABEL KONTINGENSI 3X3

Ukuran asosiasi adalah indikator statistik yang digunakan untuk mengukur kekuatan dan arah hubungan antara dua variabel, khususnya dalam konteks variabel kategorik.

7.3.1 RISK DIFFERENCE (PERBEDAAN RISIKO)

Risk Difference menunjukkan selisih risiko / probabilitas antara kelompok yang terpapar dan tidak terpapar terhadap suatu kejadian.

Rumus Risk Difference dapat dituliskan dalam formula berikut :

\[ BP = P(Y | X_1, Z) - P(Y | X_2, Z) \]

Dan kriteria uji dalam Risk Difference adalah :

Jika RD > 0, maka risiko kejadian lebih tinggi di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika RD < 0, maka risiko kejadian lebih rendah di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika RD = 0, maka tidak ada asosiasi antara dua kelompok.

7.3.2 RELATIVE RISK (RISIKO RELATIF)

Relative Risk mengukur seberapa besar risiko kejadian pada kelompok terpapar relatif terhadap kelompok tidak terpapar.

Rumus Relative Risk dapat dituliskan dalam formula berikut :

\[ RR = \frac{P(Y | X_1, Z)}{P(Y | X_2, Z)} \] Dan kriteria uji dalam Relative Risk adalah :

Jika RR > 1, maka risiko kejadian lebih tinggi di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika RR < 1, maka risiko kejadian lebih rendah di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika RR = 1, maka tidak ada asosiasi antara dua kelompok.

7.3.3 ODDS RATIO (RASIO ODDS)

Odds Ratio menunjukkan perbandingan peluang kejadian antar kelompok.

Rumus Odds Ratio dapat dituliskan dalam formula berikut :

\[ OR = \frac{P(Y | X_1, Z) / (1 - P(Y | X_1, Z))}{P(Y | X_2, Z) / (1 - P(Y | X_2, Z))} \] Dan kriteria uji dalam Odds Ratio adalah :

Jika OR > 1, maka peluang kejadian lebih tinggi di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika OR < 1, maka peluang kejadian lebih rendah di kelompok terpapar dibandingkan kelompok tidak terpapar.
Jika OR = 1, maka tidak ada asosiasi antara dua kelompok.

Contoh soal menggunakan bantuan software R studio :

Sebuah penelitian dilakukan untuk mengetahui hubungan antara jenis kelamin, usia, dan jenis penyakit (Diare/Gastritis). Penelitian ini melibatkan sejumlah pasien yang dikategorikan berdasarkan:

Jenis Kelamin: Laki-laki / Perempuan
Usia: Anak-anak / Remaja / Dewasa
Jenis Penyakit: Diare / Gastritis

Data survei yang diperoleh adalah sebagai berikut:

# Dataset contoh
data <- data.frame(
  Jenis_Kelamin = rep(c("Laki-laki", "Perempuan"), each = 3),
  Usia = rep(c("Anak-anak", "Remaja", "Dewasa"), times = 2),
  Diare = c(50,45,35,20,25,30),
  Gastritis = c(5,10,15,40,35,50)
)

data$Total_Baris <- data$Diare + data$Gastritis
total_kolom <- colSums(data[, 3:5])
total_kolom_df <- data.frame(Jenis_Kelamin = "Total Kolom", Usia = "", t(total_kolom))

data <- rbind(data, total_kolom_df)
data[3:5] <- lapply(data[3:5], as.numeric)

kable(data, caption = "Tabel Kontingensi Penyakit berdasarkan Jenis Kelamin dan Usia", 
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE) %>%
  row_spec(nrow(data), bold = TRUE, italic = FALSE)%>%
  collapse_rows(columns=1, valign="middle")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Kontingensi Penyakit berdasarkan Jenis Kelamin dan Usia
Jenis_Kelamin	Usia	Diare	Gastritis	Total_Baris
Laki-laki	Anak-anak	50	5	55
	Remaja	45	10	55
	Dewasa	35	15	50
Perempuan	Anak-anak	20	40	60
	Remaja	25	35	60
	Dewasa	30	50	80
Total Kolom		205	155	360

Risk Difference

a <- sum(data$Diare[data$Usia == "Dewasa"])
b <- sum(data$Diare[data$Usia == "Anak-anak"])
c <- sum(data$Gastritis[data$Usia == "Dewasa"])
d <- sum(data$Gastritis[data$Usia == "Anak-anak"])

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

## [1] -0.1086957

Melalui perhitungan Risk Difference tersebut, karena RD < 1, maka kelompok dewasa memiliki risiko lebih rendah terkena diare dibandingkan anak-anak.

Relative Risk

# Relative Risk
RR <- risk_dewasa / risk_anak
RR

## [1] 0.8214286

Melalui perhitungan Relative Risk tersebut, karena RR < 1, maka kelompok dewasa memiliki risiko lebih rendah terkena diare dibandingkan anak-anak.

Odds Ratio

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

## [1] 0.6428571

Melalui perhitungan Odds Ratio tersebut, karena OR < 1, maka kelompok dewasa memiliki risiko peluang lebih rendah terkena diare dibandingkan anak-anak.

7.4 INDEPENDENSI BERSYARAT

Independensi bersyarat adalah suatu keadaan dimana dua variabel acak menjadi independen jika hubungan antara keduanya menjadi tidak signifikan saat variabel ketiga dipertimbangkan.

Rumus independensi bersyarat dapat dituliskan dalam formula berikut : \[ P(X, Y | Z) = P(X | Z) \cdot P(Y | Z) \] Contoh soal menggunakan bantuan software R studio :

Di sebuah universitas, dilakukan penelitian kecil terhadap 160 mahasiswa dari dua fakultas: Fakultas Sains dan Fakultas Sosial. Penelitian ini ingin mengetahui hubungan antara kebiasaan minum kopi saat belajar (X) dan hasil ujian akhir (Y), dengan mempertimbangkan fakultas asal mahasiswa (Z). Data penelitian yang diperoleh adalah sebagai berikut:

# Dataset contoh
data <- expand.grid(
  X = c("Ya", "Tidak"),
  Y = c("Lulus", "Tidak"),
  Z = c("Sains", "Sosial"))

data$Freq <- c(25, 5, 15, 15, 20, 10, 10, 20)
table3d <- xtabs(Freq ~ X + Y + Z, data = data)

tabel_data <- data.frame(
  Kopi = c("Ya", "Ya", "Tidak", "Tidak"),
  Nilai = c("Lulus", "Tidak", "Lulus", "Tidak"),
  Sains = c(25, 15, 5, 15),
  Sosial = c(20, 10, 10, 20)
)

kable(tabel_data,
      col.names = c("Kopi", "Nilai", "Sains", "Sosial"),
      caption = "Tabel Kontingensi: Kopi, Nilai, dan Bidang Studi",
      booktabs = TRUE, escape = FALSE) %>%
  kable_styling(full_width = FALSE) %>%
  collapse_rows(columns = 1, valign = "middle")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Kontingensi: Kopi, Nilai, dan Bidang Studi
Kopi	Nilai	Sains	Sosial
Ya	Lulus	25	20
Ya	Tidak	15	10
Tidak	Lulus	5	10
Tidak	Tidak	15	20

# Chi-square
sains <- table3d[, , "Sains"]
sosial <- table3d[, , "Sosial"]

# Fakultas Sains
cat("Uji Chi-Square untuk Bidang Sains:\n")

## Uji Chi-Square untuk Bidang Sains:

uji_sains <- chisq.test(sains, correct = FALSE)
print(uji_sains)

## 
##  Pearson's Chi-squared test
## 
## data:  sains
## X-squared = 7.5, df = 1, p-value = 0.00617

# Fakultas Sosial
cat("\nUji Chi-Square untuk Bidang Sosial:\n")

## 
## Uji Chi-Square untuk Bidang Sosial:

uji_sosial <- chisq.test(sosial, correct = FALSE)
print(uji_sosial)

## 
##  Pearson's Chi-squared test
## 
## data:  sosial
## X-squared = 6.6667, df = 1, p-value = 0.009823

Melalui uji independensi bersyarat tersebut, didapatkan beberapa kesimpulan sebagai berikut :

Untuk mahasiswa fakultas sains, karena p-value bernilai 0,00617 ( p-value < 0,05) maka hipotesis nol ditolak. Dengan kata lain, kebiasaan minum kopi berpengaruh terhadap kelulusan di fakultas sains.
Untuk mahasiswa fakultas sosial, karena p-value bernilai 0,009823 ( p-value < 0,05) maka hipotesis nol ditolak. Dengan kata lain, kebiasaan minum kopi berpengaruh terhadap kelulusan di fakultas sosial.
Karena variabel kopi dan kelulusan saling berhubungan bahkan setelah memperhitungkan fakultas, maka kedua variabel tersebut tidak kondisional independen terhadap fakultas.

7.5 INFERENSI TABEL KONTINGENSI TIGA ARAH

Tabel kontingensi tiga arah digunakan untuk menganalisis hubungan antara dua variabel kategorik dengan mempertimbangkan variabel kontrol.

7.5.1 COCHRAN-MANTEL-HAENSZEL

Uji Cochran-Mantel-Haenszel (CMH) adalah metode statistik yang digunakan untuk menguji asosiasi antara dua variabel kategori dengan mengendalikan efek dari variabel ketiga yang juga bersifat kategori. Uji ini berguna ketika data disusun dalam beberapa strata atau lapisan, dan kita ingin mengetahui apakah hubungan antara dua variabel utama konsisten di seluruh strata tersebut.

Tujuan utama dari uji CMH adalah sebagai berikut :

Menentukan apakah terdapat hubungan signifikan antara dua variabel utama setelah mengontrol variabel perancu (confounding).
Menghitung odds ratio gabungan yang memperhitungkan efek dari variabel perancu.

Hipotesis dalam uji CMH adalah :

H0 = Tidak ada hubungan antara variabel kategori
H1 = Terdapat hubungan antara variabel kategori

Rumus CMH dapat dituliskan dalam formula berikut :

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

n1ik: nilai frekuensi sel baris 1 kolom 1 pada tabel parsial ke-k.

miu1ik: nilai ekspektasi sel baris 1 kolom 1 pada tabel parsial ke-k, dihitung dengan rumus: \[ \mu_{1ik} = E(n_{1ik}) = \frac{n_{1.k} \cdot n_{.1k}}{n_{..}} \]

Varians dari n1ik diberikan oleh: \[ \text{var}(n_{1ik}) = \frac{n_{1.k} \cdot n_{.2k} \cdot n_{1.k} \cdot n_{2k}}{n_{..}^2} \]

Dan kriteria uji dalam uji CMH adalah :

Jika p-value < 0,05 maka tolak H0 dan dapat disimpulkan bahwa terdapat hubungan antara variabel kategori.
Jika p-value > 0,05 maka terima H0 dan dapat disimpulkan bahwa tidak terdapat hubungan yang signifikan antara variabel kategori.

7.5.2 ODDS RATIO BERSAMA

Odds Ratio (OR) mengukur kekuatan asosiasi atau ukuran pengaruh antara dua kategori (misalnya, terkena penyakit dan tidak terkena) dalam dua kelompok atau lebih. Odds Ratio Bersama dihitung untuk menggabungkan hasil Odds Ratio dari berbagai studi atau kelompok dalam satu kesimpulan.

Rumus Odds Ratio dapat dituliskan dalam formula berikut :

\[ \hat{\theta}_{MH} = \frac{\sum_{k=1}^K \left(\frac{n_{11k} n_{22k}}{n_{\cdot k}}\right)}{\sum_{k=1}^K \left(\frac{n_{12k} n_{21k}}{n_{\cdot k}}\right)} \]

Dimana:

n11k: Frekuensi sel baris 1 kolom 1 pada tabel parsial ke-k

n12k: Frekuensi sel baris 1 kolom 2 pada tabel parsial ke-k

n21k: Frekuensi sel baris 1 pada tabel parsial ke-k

n22k: Frekuensi sel baris 2 pada tabel parsial ke-k

n..k: Total observasi dalam tabel parsial ke-k.

Standard error untuk log odds ratio dapat dituliskan dalam formula berikut :

\[ \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})(\sum n_{12k}n_{21k}/n_{..k})} + \frac{{\sum(n_{12k} + n_{212k})(n_{12k}n_{21k})}/{n_{..k}^2}}{2({\sum n_{12k}n_{21k}}/{n_{..k}})^2} \]

Interval kepercayaan untuk log odds ratio dapat dituliskan dalam formula berikut :

\[ \log(\hat{\theta}_{MH}) \pm Z_{\alpha/2}\sigma^2[\log(\hat{\theta}_{MH})] \]

7.5.3 UJI HOMOGENITAS ODDS RATIO DENGAN BRESLOW-DAY

Uji Breslow-Day digunakan untuk menguji homogenitas Odds Ratio, yaitu apakah OR yang dihitung dari beberapa kelompok atau strata adalah konsisten atau homogen. Uji ini sangat penting ketika kita ingin memastikan bahwa tidak ada perbedaan signifikan dalam ukuran efek antar kelompok atau strata yang berbeda.

Hipotesis dalam uji Breslow-Day adalah :

H0 = Tidak ada hubungan antara variabel kategori
H1 = Terdapat hubungan antara variabel kategori

Rumus Breslow-Day dapat dituliskan dalam formula berikut :

\[ X^2_{\text{HBD}} = \sum_{j=1}^{K} \frac{(a_j - \tilde{a}_j)^2}{\widehat{\mathrm{Var}}\left(a_j \mid \widehat{\mathrm{OR}}_{\text{MH}}\right)} \]

dimana:

\(a_j\) adalah jumlah kasus terpapar yang diamati dalam strata \(j\).

\(\tilde{a}_j\) adalah jumlah kasus terpapar yang dihargakan berdasarkan hipotesis nol.

\(\widehat{\mathrm{Var}}\left(a_j \mid \widehat{\mathrm{OR}}_{\mathrm{MH}}\right)\) adalah varians dari \(a_j\).

Dan kriteria uji dalam uji Breslow-Day adalah :

Jika p-value < 0,05 maka tolak H0 dan dapat disimpulkan bahwa terdapat hubungan antara variabel kategori.
Jika p-value > 0,05 maka terima H0 dan dapat disimpulkan bahwa tidak terdapat hubungan yang signifikan antara variabel kategori.

Contoh soal menggunakan bantuan software R studio :

Sebuah penelitian dilakukan untuk mengetahui hubungan antara pekerjaan, tingkat pendidikan, dan usia. Penelitian ini melibatkan sejumlah orang yang dikategorikan berdasarkan:

Status Pekerjaan: Bekerja / Tidak Bekerja
Tingkat Pendidikan: S1 / S2
Usia: 20-30 / 31-40 / 41-50

Data survei yang diperoleh adalah sebagai berikut:

# Dataset contoh
data <- expand.grid(
  Usia = c("20-30", "31-40", "41-50"),
  Status_Pekerjaan = c("Bekerja", "Tidak Bekerja"),
  Pendidikan = c("S1", "S2")
)

data$Freq <- c(40, 10, 20, 5, 50, 15, 30, 10, 35, 5, 15, 10)

tabel_data <- data.frame(
  Usia = c(rep("20-30", 2), rep("31-40", 2), rep("41-50", 2)),
  Status_Pekerjaan = rep(c("Bekerja", "Tidak Bekerja"), 3),
  S1 = c(40, 20, 50, 30, 35, 15),
  S2 = c(10, 5, 15, 10, 5, 10)
)

kable(tabel_data, 
      col.names = c("Usia", "Status Pekerjaan", "S1", "S2"), 
      caption = "Tabel Data Status Pekerjaan dan Tingkat Pendidikan Berdasarkan Usia", 
      booktabs = TRUE) %>%
  kable_styling(full_width = FALSE)

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Tabel Data Status Pekerjaan dan Tingkat Pendidikan Berdasarkan Usia
Usia	Status Pekerjaan	S1	S2
20-30	Bekerja	40	10
20-30	Tidak Bekerja	20	5
31-40	Bekerja	50	15
31-40	Tidak Bekerja	30	10
41-50	Bekerja	35	5
41-50	Tidak Bekerja	15	10

# Inferensi tabel kontingensi tiga arah
table3d <- xtabs(Freq ~ Status_Pekerjaan + Pendidikan + Usia, data = data)

or_20_30 <- (40 * 5) / (20 * 10)  
or_31_40 <- (50 * 10) / (30 * 15) 
or_41_50 <- (35 * 10) / (15 * 5)  

# Odds Ratio Bersama
OR_bersama <- (or_20_30 + or_31_40 + or_41_50) / 3
OR_bersama

## [1] 2.259259

# CMH
library(epitools)
cmh_test <- mantelhaen.test(table3d)
cmh_test

## 
##  Mantel-Haenszel chi-squared test with continuity correction
## 
## data:  table3d
## Mantel-Haenszel X-squared = 5.1437, df = 1, p-value = 0.02333
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
##  0.2536889 0.8796689
## sample estimates:
## common odds ratio 
##         0.4724005

#Breslow-Day manual
strata <- list(
  `20-30` = matrix(c(40, 20, 10, 5), nrow = 2, byrow = TRUE),
  `31-40` = matrix(c(50, 30, 15, 10), nrow = 2, byrow = TRUE),
  `41-50` = matrix(c(35, 15, 5, 10), nrow = 2, byrow = TRUE)
)
# Hitung CMH OR gabungan
mh <- mantelhaen.test(array(c(40,20,10,5, 50,30,15,10, 35,15,5,10), 
                            dim = c(2,2,3)))
or_cmh <- mh$estimate
# Fungsi untuk hitung komponen Breslow-Day
breslow_day_component <- function(mat, or_cmh) {
  a <- mat[1,1]; b <- mat[1,2]; c <- mat[2,1]; d <- mat[2,2]
  n <- a + b + c + d
  r1 <- a + b
  c1 <- a + c
  est_a <- (r1 * c1 * or_cmh) / ((or_cmh * c1) + (n - c1))
  var_a <- r1 * (n - r1) * c1 * (n - c1) * or_cmh / (( (or_cmh * c1 + (n - c1))^2 ) * (n - 1))
  Q <- (a - est_a)^2 / var_a
  return(Q)
}
# Hitung Q_BD
Q_BD <- sum(sapply(strata, breslow_day_component, or_cmh = or_cmh))
Q_BD

## [1] 29.14263

# Bandingkan dengan chi-square
p_val <- 1 - pchisq(Q_BD, df = 2)
p_val

## [1] 4.69633e-07

# Hasil 
library(knitr)
library(kableExtra)

# Simpan hasil CMH test
cmh_pval <- cmh_test$p.value
cmh_stat <- cmh_test$statistic
cmh_df   <- cmh_test$parameter

# Buat tabel hasil inferensi
hasil_uji <- data.frame(
  Uji = c("Breslow-Day", "Cochran-Mantel-Haenszel", "Rata-rata Odds Ratio"),
  Statistik = c(round(Q_BD, 3), round(cmh_stat, 3), round(OR_bersama, 3)),
  `Derajat Bebas` = c(2, cmh_df, NA),
  `p-value` = c(format.pval(1 - pchisq(Q_BD, df = 2), digits = 3),
                format.pval(cmh_pval, digits = 3),
                NA)
)

# Tampilkan tabel
kable(hasil_uji,
      caption = "Hasil Uji Breslow-Day, CMH, dan Odds Ratio Bersama",
      booktabs = TRUE,
      col.names = c("Uji", "Statistik", "Derajat Bebas", "p-value")) %>%
  kable_styling(full_width = FALSE, position = "center")

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")
## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Hasil Uji Breslow-Day, CMH, dan Odds Ratio Bersama
Uji	Statistik	Derajat Bebas	p-value
Breslow-Day	29.143	2	4.7e-07
Cochran-Mantel-Haenszel	5.144	1	0.0233
Rata-rata Odds Ratio	2.259	NA	NA

Melalui inferensi tabel kontingensi 3 arah tersebut, didapatkan beberapa kesimpulan sebagai berikut :

Untuk uji CMH, didapatkan nilai p-value senilai 0,0233 (p-value < 0,05) sehingga hipotesis nol ditolak. Dengan kata lain, terdapat asosiasi antara status pekerjaan dan tingkat pendidikan setelah mengendalikan perbedaan usia.
Untuk Odds Ratio gabungan, diperoleh nilai OR sebesar 2,259 (OR > 1) sehingga dapat disimpulkan bahwa orang dengan pendidikan S1 memiliki 2,26 kali peluang bekerja dibandingkan S2, tanpa memperhitungkan variabilitas antar usia.
Untuk uji Breslow-Day, diperoleh nilai p-value sebesar 4,7x10^-7 (p-value < 0,05) sehingga dapat disimpulkan bahwa odds ratio antar kelompok usia tidak homogen. Selain itu, efek pendidikan terhadap pekerjaan tidak homogen pula di setiap strata.

8 GENERALIZED LINEAR MODEL (GLM)

Generalized Linear Model (GLM) merupakan perluasan dari model regresi linear klasik. GLM memungkinkan kita untuk memodelkan data di mana:

Variabel respons tidak berdistribusi normal.
Hubungan antara variabel prediktor dan ekspektasi dari respons tidak linear.

GLM terdiri dari tiga komponen utama:

Distribusi dari exponential family untuk variabel respons.
Fungsi link yang menghubungkan ekspektasi dari variabel respons ke kombinasi linear dari prediktor.
Fungsi linear prediktor:
\[ \eta = X \beta \] ## EXPONENTIAL FAMILY Distribusi termasuk dalam exponential family jika dapat ditulis dalam bentuk:

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

Contoh distribusi yang termasuk dalam exponential family:

Distribusi Normal
Distribusi Binomial
Distribusi Poisson
Distribusi Gamma

Contoh Pembuktian Distribusi Binomial

Fungsi probabilitas distribusi binomial dapat diformulasikan dalam rumus berikut:

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

Apabila ditulis ulang ke bentuk exponential family:

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

Dengan substitusi:

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

Maka distribusi Binomial termasuk dalam exponential family.

8.1 MODEL REGRESI LOGISTIK

Persamaan regresi logistik menyerupai regresi linear, di mana nilai input dikombinasikan secara linear dengan koefisien (bobot) untuk menghasilkan prediksi. Namun, regresi logistik membatasi hasil prediksi menjadi nilai biner, yaitu 0 atau 1, dengan menggunakan fungsi aktivasi sigmoid.

Output yang diprediksi oleh model terletak dalam rentang antara 0 hingga 1 dan mengikuti bentuk kurva S (S-shaped).

Regresi logistik menganalisis hubungan antara satu atau lebih variabel independen dan mengklasifikasikan data ke dalam kelas-kelas diskrit. Model ini banyak digunakan dalam pemodelan prediktif, di mana model memperkirakan probabilitas matematis apakah suatu entitas termasuk ke dalam kategori tertentu atau tidak.

Sebagai contoh, angka 0 dapat mewakili kelas negatif, dan angka 1 mewakili kelas positif. Regresi logistik biasanya digunakan untuk masalah klasifikasi biner, di mana variabel hasil hanya memiliki dua kemungkinan kategori (0 dan 1).

Contoh penerapan klasifikasi biner yaitu dalam e-commerce, pelanggan yang memberi rating rendah (1) atau tinggi (0) dikaitkan dengan variabel seperti waktu pengiriman, kualitas produk, dan kemudahan pengembalian. Model logistik memprediksi probabilitas pelanggan memberi rating buruk. Bisa digunakan untuk sistem rekomendasi logistik atau QC (quality control).

Persamaan dan asumsi regresi logistik adalah sebagai berikut : Regresi logistik menggunakan fungsi logistik (sigmoid) untuk mengubah output prediksi menjadi probabilitas antara 0 dan 1. Model ini sangat sesuai untuk klasifikasi biner.

Fungsi Sigmoid :

\[ f(x) = \frac{1}{1 + e^{-x}} \]

Jika: \(f(x) > 0.5 \Rightarrow \text{Kelas 1 (positif)}\)
\(f(x) < 0.5 \Rightarrow \text{Kelas 0 (negatif)}\)

Contoh Simulasi dengan data buatan dan bantuan software R Studio :

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

# Buat data frame
data <- data.frame(x = x, y = y, prob = prob)

# Visualisasi menggunakan base R
plot(x, y, pch = 16, col = "gray60",
     xlab = "X", ylab = "Y / Probabilitas",
     main = "Simulasi Regresi Logistik dengan Kurva Sigmoid")

lines(x, prob, col = "blue", lwd = 2)
abline(h = 0.5, col = "red", lty = 2)

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

Kurva sigmoid dalam regresi logistik menunjukkan hubungan non-linear antara variabel prediktor dan probabilitas output. Model ini cocok digunakan ketika target yang diprediksi adalah variabel biner, misalnya ya/tidak, sukses/gagal, positif/negatif, dll.).

Spesifikasi model regresi logistik : - Fungsi link (logit): \((\mu) = \log\left( \frac{\mu}{1 - \mu} \right)\)

Model regresi logistik: \(\log\left( \frac{\mu}{1 - \mu} \right) = X\beta\)
Fungsi inverse link: \(\mu = \frac{\exp(X\beta)}{1 + \exp(X\beta)}\)

8.1.1 ESTIMASI PARAMETER

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

Log-likelihood fungsi untuk regresi logistik: () = _{i=1}^{n} dimana : \(\pi_i = \frac{\exp(x_i^\top \beta)}{1 + \exp(x_i^\top \beta)}\)

Contoh Simulasi dengan data buatan dan bantuan software R Studio :

1. ESTIMASI DENGAN ITERASI NEWTON-RAPHSON / ALGORITMA FISHER SCORING

set.seed(123)
n <- 150
x <- rnorm(n)
log_odds <- -0.3 + 2 * x
prob <- 1 / (1 + exp(-log_odds))
y <- rbinom(n, 1, prob)
data <- data.frame(y = y, x = x)

# Fit model regresi logistik
model <- glm(y ~ x, data = data, family = binomial)
summary(model)

## 
## Call:
## glm(formula = y ~ x, family = binomial, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -0.4437     0.2024  -2.192   0.0284 *  
## x             1.6742     0.2940   5.694 1.24e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 203.41  on 149  degrees of freedom
## Residual deviance: 150.08  on 148  degrees of freedom
## AIC: 154.08
## 
## Number of Fisher Scoring iterations: 5

# Visualisasi
plot(data$x, data$y, pch = 16, col = "grey")
curve(predict(model, newdata = data.frame(x = x), type = "response")[order(x)],
add = TRUE, col = "blue", lwd = 2)

SIMULASI DATA BINER

set.seed(2025)
n <- 160
x <- rnorm(n, mean = 0.5)
log_odds <- -0.6 + 1.8 * x
prob <- 1 / (1 + exp(-log_odds))
y <- rbinom(n, 1, prob)
data <- data.frame(x = x, y = y)

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

## 
## Call:
## glm(formula = y ~ x, family = binomial, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -0.4020     0.2125  -1.892   0.0585 .  
## x             1.6791     0.2744   6.119 9.43e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 220.20  on 159  degrees of freedom
## Residual deviance: 158.87  on 158  degrees of freedom
## AIC: 162.87
## 
## Number of Fisher Scoring iterations: 5

exp(coef(model))

## (Intercept)           x 
##   0.6689575   5.3609443

library(ggplot2)
data$pred <- predict(model, type = "response")

ggplot(data, aes(x = x, y = y)) +
  geom_point(alpha = 0.4, color = "gray40") +
  geom_line(aes(y = pred), color = "blue", size = 1.2) +
  labs(title = "Kurva Logit dari Model Regresi Logistik",
       x = "X (Prediktor)", y = "Probabilitas Prediksi") +
  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.

data$pred_class <- ifelse(data$pred > 0.5, 1, 0)
table(Prediksi = data$pred_class, Aktual = data$y)

##         Aktual
## Prediksi  0  1
##        0 50 19
##        1 22 69

Beberapa kesimpulan yang diperoleh melalui uji regresi logistik di atas adalah : - Koefisien untuk Intercept adalah -0,4020, yang menunjukkan bahwa ketika variabel x=0, log odds untuk y=1 adalah -0,4020. Koefisien untuk x adalah 1,6791, yang berarti setiap kenaikan satu unit pada x meningkatkan log odds untuk y=1 sebesar 1,6791. Ini menunjukkan bahwa x berpengaruh signifikan terhadap probabilitas y=1.

Nilai p untuk x adalah 9,43 x 10^-10 yang menunjukkan bahwa x memiliki pengaruh yang sangat kuat terhadap probabilitas y=1.
Berdasarkan confusion matrix, model ini memiliki akurasi sekitar 74%, dengan kemampuan yang cukup baik dalam mendeteksi kelas positif (kelas 1).

8.2 MODEL REGRESI POISSON

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

Distribusi Poisson memiliki fungsi probabilitas:

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

Bentuk exponential family:

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

Dengan:

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

Fungsi link kanonik untuk distribusi Poisson adalah fungsi logaritma:

\[ g(\mu) = \log(\mu) \Rightarrow \log(\mu_i) = \mathbf{x}_i^\top \boldsymbol{\beta} \]

Inverse link-nya:

\[ \mu_i = \exp(\mathbf{x}_i^\top \boldsymbol{\beta}) \]

8.2.1 ESTIMASI PARAMETER

Estimasi parameter \(\boldsymbol{\beta}\) dilakukan dengan metode Maximum Likelihood Estimation (MLE), dengan fungsi log-likelihood:

\[ \ell(\beta) = \sum_{i=1}^{n} [y_i \mathbf{x}_i^\top \beta - \exp(\mathbf{x}_i^\top \beta) - \log(y_i!)] \]

Contoh Simulasi dengan data buatan dan bantuan software R Studio :

set.seed(100)
n <- 150
x <- rnorm(n)
lambda <- exp(0.2 + 0.5 * x)
y <- rpois(n, lambda)
data <- data.frame(y = y, x = x)

poisson_model <- glm(y ~ x, data = data, family = poisson)
summary(poisson_model)

## 
## Call:
## glm(formula = y ~ x, family = poisson, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.22776    0.07591   3.001  0.00269 ** 
## x            0.51029    0.07021   7.268 3.64e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 202.45  on 149  degrees of freedom
## Residual deviance: 150.93  on 148  degrees of freedom
## AIC: 424.15
## 
## Number of Fisher Scoring iterations: 5

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

# Diagnostik Overdispersi
dispersion <- sum(residuals(poisson_model, type = "pearson")^2) / poisson_model$df.residual
dispersion

## [1] 0.9373263

Beberapa kesimpulan yang diperoleh melalui uji regresi logistik di atas adalah : - Koefisien untuk Intercept adalah 0,227, yang menunjukkan bahwa ketika variabel x=0, log odds untuk y=1 adalah 0,227. Koefisien untuk x adalah 0,51, yang berarti setiap kenaikan satu unit pada x meningkatkan log odds untuk y=1 sebesar 0,51. Ini menunjukkan bahwa x berpengaruh signifikan terhadap probabilitas y=1.

Nilai p untuk x adalah 3,64 x 10^-13 yang menunjukkan bahwa x memiliki pengaruh yang sangat kuat terhadap probabilitas y=1.
Berdasarkan nilai diagnostik overdispersi, karena nilai < 1, model ini tidak mengalami gejala overdispersi.

# Dataset
data("warpbreaks")
warpbreaks_mod <- glm(breaks ~ wool + tension, data = warpbreaks, family = poisson)
summary(warpbreaks_mod)

## 
## Call:
## glm(formula = breaks ~ wool + tension, family = poisson, data = warpbreaks)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  3.69196    0.04541  81.302  < 2e-16 ***
## woolB       -0.20599    0.05157  -3.994 6.49e-05 ***
## tensionM    -0.32132    0.06027  -5.332 9.73e-08 ***
## tensionH    -0.51849    0.06396  -8.107 5.21e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 297.37  on 53  degrees of freedom
## Residual deviance: 210.39  on 50  degrees of freedom
## AIC: 493.06
## 
## Number of Fisher Scoring iterations: 4

exp(coef(warpbreaks_mod))

## (Intercept)       woolB    tensionM    tensionH 
##  40.1235380   0.8138425   0.7251908   0.5954198

warpbreaks$predicted <- predict(warpbreaks_mod, type = "response")

library(ggplot2)
ggplot(warpbreaks, aes(x = tension, y = breaks, color = wool)) +
  geom_jitter(width = 0.2, alpha = 0.6) +
  geom_point(aes(y = predicted), shape = 18, size = 3, color = "black") +
  facet_wrap(~wool) +
  labs(title = "Prediksi Patahan Benang Berdasarkan Wol dan Ketegangan",
       x = "Tingkat Ketegangan", y = "Jumlah Patahan", color = "Jenis Wol") +
  theme_minimal()

# Plot residuals
plot(poisson_model$residuals, main ="Residual Plot", ylab = "Residual", xlab = "Index", pch = 19, col = "steelblue")
abline(h = 0, col = "red", lty = 2)

Beberapa kesimpulan yang diperoleh melalui uji regresi logistik di atas adalah : - Koefisien untuk Intercept adalah 3,69, yang menunjukkan bahwa ketika variabel x=0, log odds untuk y=1 adalah 3,69. Koefisien untuk x menunjukkan nilai < 1. Ini menunjukkan bahwa x berpengaruh signifikan terhadap probabilitas y=1 dengan arah berkelbalikan, yaitu semakin x bertambah maka semakin berkurang jumlah variabel y.

9 INFERENSI GENERALIZED LINEAR MODEL (GLM)

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

9.1 EKSPEKTASI DAN VARIANS

Untuk distribusi exponential family, ekspektasi dan varians dari Yi adalah:

\[ \mathbb{E}[Y_i] = \mu_i,\quad \mathrm{Var}(Y_i) = \phi V(\mu_i) \]

Contoh fungsi varians:

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

Dengan asumsi distribusi asimptotik:

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

# Simulasi regresi Poisson
set.seed(321)
x <- rnorm(120)
mu <- exp(0.4 + 0.6 * x)
y <- rpois(120, mu)
model <- glm(y ~ x, family = poisson)

summary(model)

## 
## Call:
## glm(formula = y ~ x, family = poisson)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.35353    0.08213   4.305 1.67e-05 ***
## x            0.68246    0.07639   8.934  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 220.88  on 119  degrees of freedom
## Residual deviance: 139.70  on 118  degrees of freedom
## AIC: 374.52
## 
## Number of Fisher Scoring iterations: 5

9.1.1 ESTIMASI PARAMETER

Pendekatan numerik seperti Newton-Raphson dan Fisher Scoring digunakan dua metode berikut untuk mempercepat konvergensi estimasi, yakni:

Newton-Raphson \[ \beta^{(t+1)} = \beta^{(t)} - \left[ H(\beta^{(t)}) \right]^{-1} U(\beta^{(t)}) \]

dan IRLS (Iteratively Reweighted Least Squares).

9.2 DIAGNOSTIK MODEL GLM

Statistik yang digunakan untuk menilai kesesuaian model adalah sebagai berikut : - Deviance \[ D = 2 \sum_{i=1}^{n} \left[ y_i \log\left(\frac{y_i}{\hat{\mu}_i}\right) - (y_i - \hat{\mu}_i) \right] \]

Pearson Chi-square \[ X^2 = \sum_{i=1}^{n} \frac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i} \]

deviance(model)

## [1] 139.7009

sum(residuals(model, type = "pearson")^2)

## [1] 134.4544

# Plot Residual
plot(residuals(model), main = "Residual Plot", ylab = "Residual", col = "blue", pch = 19)
abline(h = 0, col = "red", lty = 2)

9.3 DETAIL METODE ESTIMASI DAN INFERENSI REGRESI LOGISTIK

Fungsi model logistik : \[ \pi(x) = \frac{\exp(\beta_0 + \beta_1 x)}{1 + \exp(\beta_0 + \beta_1 x)} \] Fungsi log-likelihood yang akan dimaksimumkan adalah:

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

9.3.1 ESTIMASI DENGAN NEWTON-RAPHSON

Model probabilitas dapat dituliskan sebagai:

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

Fungsi log-likelihood:

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

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

X <- as.matrix(X)
beta <- rep(0, ncol(X))

tol <- 1e-6
max_iter <- 100
for (i in 1:max_iter) {
eta <- X %*% beta
pi_hat <- 1 / (1 + exp(-eta))
W <- diag(as.numeric(pi_hat * (1 - pi_hat)))
z <- eta + solve(W) %*% (y - pi_hat)
beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z
if (sum(abs(beta_new - beta)) < tol) break
beta <- beta_new
}
beta

##        [,1]
##   -1.516679
## x  2.155658

9.4 UJI WALD

Uji Wald bertujuan untuk menguji signifikansi parameter \(\beta_{j}\) dalam model regresi logistik.

Hipotesis untuk uji Wald adalah : H0 = \(\beta_{j}\) = 0 H1 = \(\beta_{j}\) ≠ 0

Dengan rumus statistik Wald : \[ W = \left( \frac{\hat{\beta}_j}{\mathrm{SE}(\hat{\beta}_j)} \right)^2 \sim \chi^2_1 \]

# Dataset 
set.seed(123)
n <- 100
x <- rnorm(n)
log_odds <- -0.5 + 1.2 * x
p <- 1 / (1 + exp(-log_odds))
y <- rbinom(n, 1, p)
data <- data.frame(x, y)
model <- glm(y ~ x, data = data, family = binomial)
summary(model)

## 
## Call:
## glm(formula = y ~ x, family = binomial, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -0.3097     0.2296  -1.349    0.177    
## x             1.2663     0.3080   4.111 3.94e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 137.99  on 99  degrees of freedom
## Residual deviance: 114.76  on 98  degrees of freedom
## AIC: 118.76
## 
## Number of Fisher Scoring iterations: 4

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

Z <- beta_hat / se_beta
Z

##        x 
## 4.110965

Wald_stat <- Z^2
Wald_stat

##        x 
## 16.90003

p_value <- 1 - pchisq(Wald_stat, df = 1)
p_value

##            x 
## 3.940095e-05

Berdasarkan uji Wald tersebut, karena p-value bernilai 3,94 x 10^-5 (p-value < 0,05) maka H0 ditolak. Dengan kata lain, variabel x berpengaruh secara signifikan.

9.5 UJI LIKELIHOOD RATIO

Rumus statistik likelihood ratio : \[ G^2 = -2 \left( \ell_{\text{null}} - \ell_{\text{full}} \right) \sim \chi^2_{p - q} \]

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

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

## Analysis of Deviance Table
## 
## Model 1: y ~ 1
## Model 2: y ~ x
##   Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
## 1        99     137.99                          
## 2        98     114.76  1   23.229 1.438e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

9.6 AIC DAN BIC

Semakin kecil AIC, semakin baik model.

AIC(model)

## [1] 118.7598

Alternatif AIC, mengukur kompleksitas model.

BIC(model)

## [1] 123.9701

9.7 IRLS

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

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

# Inisialisasi
beta <- c(0, 0)
tol <- 1e-6
max_iter <- 100
for (i in 1:max_iter) {
eta <- X %*% beta
lambda <- exp(eta)
W <- diag(as.numeric(lambda))
z <- eta + (y - lambda) / lambda
beta_new <- solve(t(X) %*% W %*% X) %*% t(X) %*% W %*% z
if (sum(abs(beta_new - beta)) < tol) {
cat("Konvergen pada iterasi ke-", i, "\n")
break
}
beta <- beta_new
}

## Konvergen pada iterasi ke- 8

beta

##        [,1]
##   0.4494951
## x 0.8600048

# Perbandingan dengan GLM
model_glm <- glm(y ~ x, family = poisson)
summary(model_glm)

## 
## Call:
## glm(formula = y ~ x, family = poisson)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.44950    0.08872   5.066 4.05e-07 ***
## x            0.86000    0.07463  11.523  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 245.05  on 99  degrees of freedom
## Residual deviance: 106.78  on 98  degrees of freedom
## AIC: 325.76
## 
## Number of Fisher Scoring iterations: 5

Berdasarkan hasil uji di atas, hasil manual IRLS sangat mendekari hasil glm() di R studio**, yakni koefisien beta yang bernilai 0,44 dan 0,86.

10 REGRESI LOGISTIK DENGAN PREDIKTOR NOMINAL, ORDINAL, DAN RASIO

Regresi logistik merupakan teknik statistik yang digunakan untuk menjelaskan keterkaitan antara variabel respons dikotomis (dua kemungkinan hasil) dengan satu atau lebih variabel prediktor. Variabel-variabel prediktor yang digunakan dalam model ini bisa berasal dari berbagai skala pengukuran:

Nominal: Merupakan variabel kategorik tanpa urutan yang inheren antar kategorinya. Misalnya, jenis kelamin (laki-laki atau perempuan). Untuk keperluan analisis regresi logistik, kategori ini biasanya diubah menjadi variabel dummy, contohnya dengan mengkodekan perempuan sebagai 1 dan laki-laki sebagai 0.
Ordinal: Variabel ini memiliki urutan logis antar kategori, tetapi interval antara kategori tidak harus sama. Contohnya adalah jenjang pendidikan (SMA, Sarjana, Master, Doktor). Dalam regresi logistik, variabel ordinal dapat diolah seperti data nominal (menggunakan dummy untuk tiap level), atau diperlakukan seperti data numerik dengan memberikan skor angka untuk mencerminkan urutan, seperti SMA = 1, Sarjana = 2, dan seterusnya.
Rasio: Merupakan variabel numerik kontinu yang memiliki nol mutlak dan rasio antar nilai yang bermakna, seperti pendapatan dalam satuan mata uang.

10.1 SIMULASI DATA

Sebuah startup pertanian digital ingin mengetahui faktor-faktor yang memengaruhi keberhasilan petani dalam meningkatkan hasil panen setelah menggunakan aplikasi pertanian berbasis AI selama 6 bulan. Data dikumpulkan dari 1000 petani di berbagai daerah, dan mencakup informasi berikut:

library(tibble)

# Set seed untuk replikasi
set.seed(39)

# Jumlah petani
n <- 1000

# Variabel 1: Jenis lahan (nominal)
land_type <- sample(c("Sawah", "Ladang"), n, replace = TRUE)

# Variabel 2: Skala usaha tani (ordinal)
farm_scale <- sample(c("Kecil", "Sedang", "Besar", "Sangat Besar"), 
                     n, replace = TRUE, 
                     prob = c(0.4, 0.3, 0.2, 0.1))

# Variabel 3: Jumlah pelatihan online (rasio, bilangan positif)
training_count <- rpois(n, lambda = 5)  # Poisson agar integer, rata-rata 5 pelatihan

# Model logit (koefisien disesuaikan)
logit_p <- -1 + 
  0.6 * (land_type == "Sawah") + 
  0.7 * as.numeric(factor(farm_scale, levels = c("Kecil", "Sedang", "Besar", "Sangat Besar"), ordered = TRUE)) + 
  0.2 * training_count

# Probabilitas keberhasilan panen
p <- 1 / (1 + exp(-logit_p))

# Variabel respons: keberhasilan panen (0 = tidak meningkat, 1 = meningkat)
success <- rbinom(n, size = 1, prob = p)

# Gabungkan ke dalam dataset
sim_farm <- tibble(success, land_type, farm_scale, training_count)

# Lihat 6 baris pertama
head(sim_farm)

## # A tibble: 6 × 4
##   success land_type farm_scale   training_count
##     <int> <chr>     <chr>                 <int>
## 1       1 Sawah     Kecil                     8
## 2       1 Sawah     Kecil                     8
## 3       1 Ladang    Sangat Besar              4
## 4       1 Ladang    Kecil                     5
## 5       1 Ladang    Besar                     2
## 6       1 Sawah     Besar                     9

10.2 EKSPLORASI DATA

sim_farm %>%
dplyr::group_by(success) %>%
dplyr::summarise(
  Jumlah = dplyr::n(),
Jumlah_pelatihan = mean(training_count)
)

## # A tibble: 2 × 3
##   success Jumlah Jumlah_pelatihan
##     <int>  <int>            <dbl>
## 1       0    179             4.54
## 2       1    821             5.15

10.3 PERLAKUAN VARIABEL ORDINAL

10.3.1 TREAT SEBAGAI NOMINAL (DUMMY)

library(dplyr)
sim_farm_nominal <- sim_farm %>%
mutate(
farm_scale = factor(farm_scale, levels = c("Kecil", "Sedang", "Besar", "Sangat Besar"))
)
model_nominal <- glm(success ~ land_type + farm_scale + training_count, data = sim_farm_nominal, family = binomial)
summary(model_nominal)

## 
## Call:
## glm(formula = success ~ land_type + farm_scale + training_count, 
##     family = binomial, data = sim_farm_nominal)
## 
## Coefficients:
##                        Estimate Std. Error z value Pr(>|z|)    
## (Intercept)             0.38450    0.23473   1.638 0.101409    
## land_typeSawah          0.42033    0.17082   2.461 0.013870 *  
## farm_scaleSedang        0.09951    0.18538   0.537 0.591399    
## farm_scaleBesar         1.16801    0.28962   4.033 5.51e-05 ***
## farm_scaleSangat Besar  1.33803    0.41301   3.240 0.001196 ** 
## training_count          0.13422    0.03912   3.431 0.000601 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 939.75  on 999  degrees of freedom
## Residual deviance: 891.85  on 994  degrees of freedom
## AIC: 903.85
## 
## Number of Fisher Scoring iterations: 5

Interpretasi Koefisien

Intersep (+0,3845)
- Ini adalah log-odds dasar untuk petani berlahan Ladang, dengan skala usaha Kecil, dan jumlah pelatihan = 0.
- Tidak signifikan (p = 0,101409), artinya baseline ini tidak berbeda secara signifikan dari peluang sukses 50%.
- Odds Ratio = exp(0.3845) ≈ 1,47: Peluang sukses 47% lebih tinggi dari baseline log-odds nol, tapi tidak signifikan secara statistik.
Sawah (+0,42033)
- Petani yang bertani di lahan Sawah memiliki log-odds sukses lebih tinggi 0,42 dibanding petani di lahan Ladang.
- Signifikan di level 5% (p = 0,01387), menunjukkan bahwa jenis lahan Sawah secara statistik meningkatkan peluang sukses.
- Odds Ratio = exp(0,42033) ≈ 1,52: Peluang sukses petani sawah sekitar 52% lebih tinggi dari petani ladang.
Sedang (+0,09951)
- Petani dengan skala usaha Sedang hanya memiliki log-odds sedikit lebih tinggi dari petani dengan skala Kecil, namun tidak signifikan (p = 0,591399).
- Odds Ratio = exp(0,09951) ≈ 1,10: Kenaikan peluang sekitar 10%, tidak signifikan.
Besar (+1,16801)
- Petani dengan skala usaha Besar memiliki kenaikan log-odds sebesar 1,17 dibanding skala Kecil.
- Sangat signifikan (p < 0,001), artinya skala usaha Besar sangat berkontribusi terhadap keberhasilan panen.
- Odds Ratio = exp(1,16801) ≈ 3,22: Peluang sukses lebih dari tiga kali lipat dibanding petani skala kecil.
Sangat Besar (+1,33803)
- Petani dengan skala usaha Sangat Besar juga memiliki log-odds lebih tinggi dari petani skala Kecil.
- Sangat signifikan (p = 0,001196).
- Odds Ratio = exp(1,33803) ≈ 3,81: Peluang sukses meningkat hampir 4 kali lipat dibanding petani skala kecil.
Jumlah Pelatihan (+0,13422)
- Setiap tambahan satu sesi pelatihan online meningkatkan log-odds sukses sebesar 0,134.
- Sangat signifikan (p = 0,000601).
- Odds Ratio = exp(0,13422) ≈ 1,14: Artinya, setiap sesi pelatihan tambahan meningkatkan peluang sukses sekitar 14%.
Goodness-of-Fit
- Null deviance = 939,75: deviance dari model tanpa prediktor.
- Residual deviance = 891,85: deviance dari model yang sudah mempertimbangkan prediktor.
```
-    Penurunan deviance menunjukkan bahwa model dengan prediktor **lebih informatif** daripada model tanpa prediktor.
```
- AIC = 903,85
Signifikansi Model
- Jenis lahan (Sawah), skala usaha (Besar dan Sangat Besar), serta jumlah pelatihan merupakan prediktor signifikan terhadap keberhasilan panen.
- Skala usaha Sedang tidak signifikan, kemungkinan karena terlalu dekat dengan baseline (Kecil).
- Model secara keseluruhan menunjukkan hubungan yang bermakna antara prediktor dan peluang sukses panen.
Kesimpulan Praktis
- Petani dengan lahan Sawah, skala usaha yang lebih besar, dan yang mengikuti lebih banyak pelatihan memiliki peluang lebih tinggi untuk sukses setelah menggunakan aplikasi pertanian digital.
- Program pelatihan berdampak nyata, sehingga startup sebaiknya memperluas jangkauan pelatihan ini.

10.3.2 TREAT SEBAGAI RASIO (NUMERIC RANK)

sim_farm_numeric <- sim_farm %>%
mutate(
farm_scale_numeric = case_when(
farm_scale == "Kecil" ~ 1,
farm_scale == "Sedang" ~ 2,
farm_scale == "Besar" ~ 3,
farm_scale == "Sangat Besar" ~ 4
)
)
model_numeric <- glm(success ~ land_type + farm_scale_numeric + training_count, data = sim_farm_numeric, family = binomial)
summary(model_numeric)

## 
## Call:
## glm(formula = success ~ land_type + farm_scale_numeric + training_count, 
##     family = binomial, data = sim_farm_numeric)
## 
## Coefficients:
##                    Estimate Std. Error z value Pr(>|z|)    
## (Intercept)        -0.15901    0.27872  -0.571 0.568328    
## land_typeSawah      0.44028    0.17028   2.586 0.009723 ** 
## farm_scale_numeric  0.45824    0.09707   4.721 2.35e-06 ***
## training_count      0.13074    0.03891   3.360 0.000779 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 939.75  on 999  degrees of freedom
## Residual deviance: 897.76  on 996  degrees of freedom
## AIC: 905.76
## 
## Number of Fisher Scoring iterations: 4

Interpretasi Koefisien

Intersep (-0,15901)
- Ini adalah log-odds dasar untuk petani yang lahan pertaniannya Ladang, dengan skala usaha = 1 (kecil), dan jumlah pelatihan = 0.
- Tidak signifikan (p = 0,568328), artinya baseline ini tidak berbeda secara signifikan dari probabilitas 50%.
- Odds Ratio = exp(–0,15901) ≈ 0,85: Peluang sukses dasar sedikit lebih kecil, tetapi tidak signifikan.
Sawah (+0,44028)
- Petani yang memiliki lahan bertani di sawah memiliki log-odds keberhasilan lebih tinggi 0,44 dibandingkan petani yang lahannya bukan sawah.
- p = 0,009723 (signifikan), menunjukkan jenis lahan sawah berasosiasi dengan peluang keberhasilan yang lebih tinggi.
- Odds Ratio = exp(0,44028) ≈ 1,55 → Peluang keberhasilan bertani di sawah sekitar 55% lebih tinggi daripada bukan sawah.
farm_scale_numeric (+0,45824)
- Setiap peningkatan satu unit skala usaha (misalnya dari kecil ke sedang, atau sedang ke besar) meningkatkan log-odds keberhasilan sebesar 0,458.
- p = 2,35x(10^-6) (sangat signifikan), artinya semakin besar skala usaha, semakin besar pula peluang keberhasilan panen.
- Odds Ratio = exp(0,45824) ≈ 1,58 → Setiap kenaikan satu tingkat skala usaha meningkatkan peluang keberhasilan sekitar 58%.
training_count (+0,13074)
- Setiap tambahan 1 kali pelatihan yang diikuti oleh petani, meningkatkan log-odds keberhasilan sebesar 0,131.
- p = 0,000779 (sangat signifikan), menunjukkan bahwa pelatihan yang lebih banyak memberikan kontribusi positif terhadap keberhasilan.
- Odds Ratio = exp(0.13074) ≈ 1,14. Setiap pelatihan tambahan meningkatkan peluang keberhasilan sebesar 14%.
Goodness-of-Fit
- Null deviance = 939,75: deviance dari model tanpa prediktor.
- Residual deviance = 897,76: deviance dari model yang sudah mempertimbangkan prediktor.
- Penurunan deviance menunjukkan bahwa model dengan prediktor lebih informatif daripada model tanpa prediktor.
- AIC = 905,76
Signifikansi Model
- Jenis lahan (Sawah), skala usaha (Besar dan Sangat Besar), serta jumlah pelatihan merupakan prediktor signifikan terhadap keberhasilan panen.
- Model secara keseluruhan menunjukkan hubungan yang bermakna antara prediktor dan peluang sukses panen.
Kesimpulan Praktis
- Bertani di sawah, memiliki skala usaha yang lebih besar, dan mengikuti lebih banyak pelatihan berkontribusi signifikan terhadap keberhasilan panen.

10.3.3 PERBANDINGAN MODEL

list(
AIC_Nominal = AIC(model_nominal),
AIC_Numeric = AIC(model_numeric)
)

## $AIC_Nominal
## [1] 903.8474
## 
## $AIC_Numeric
## [1] 905.7647

Interpretasi : Model AIC kecil lebih baik, sehingga treat sebagai nominal lebih dianjurkan.

10.3.4 GOODNESS-OF-FIT MODEL

nullmod <- glm(success ~ 1, data = sim_farm, 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.05097127 (df=6)
## 
## $McFadden_R2_Numeric
## 'log Lik.' 0.0446746 (df=4)

Kesimpulan :

Model ordinal sebagai nominal menjelaskan sekitar 5.1% variasi log-likelihood dibanding model tanpa prediktor. Sementara, model ordinal sebagai sedikit lebih buruk dibanding model pertama karena hanya menjelaskan sekitar 4,6% variasi log-likelihood.
Dengan kata lain, perlakuan ordinal sebagai numerik mungkin kurang fleksibel untuk menangkap variasi antar kategori skala usaha tani.

10.3.5 VISUALISASI PREDIKSI

sim_farm_nominal <- sim_farm_nominal %>% mutate(predicted = predict(model_nominal, type = "response"))
sim_farm_numeric <- sim_farm_numeric %>% mutate(predicted = predict(model_numeric, type = "response"))

# Plot untuk model nominal
library(ggplot2)
library(dplyr)
sim_farm_nominal %>%
  ggplot(aes(x = training_count, y = predicted, color = farm_scale)) +
  geom_point(alpha = 0.6) +
  labs(
    title = "Prediksi Probabilitas (Ordinal sebagai Nominal)",
    x = "Pelatihan",
    y = "Prediksi Probability"
  ) +
  theme_minimal()

# Plot untuk model numeric
sim_farm_numeric %>%
  ggplot(aes(x = training_count, y = predicted, color = as.factor(farm_scale))) +
  geom_point(alpha = 0.6) +
  labs(
    title = "Prediksi Probabilitas (Ordinal sebagai Numerik)",
    x = "Pelatihan",
    y = "Prediksi Probability"
  ) +
  theme_minimal()

10.3.6 RINGKASAN KOEFISIEN MODEL

library(knitr)
library(kableExtra)
library(broom)

## Warning: package 'broom' was built under R version 4.4.3

# Ringkasan model nominal
tidy(model_nominal) %>%
  kable(format = "html", booktabs = TRUE, caption = "Ringkasan Koefisien Model dengan Ordinal sebagai Nominal") %>%
  kable_styling(latex_options = c("hold_position", "striped"))

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Ringkasan Koefisien Model dengan Ordinal sebagai Nominal
term	estimate	std.error	statistic	p.value
(Intercept)	0.3844967	0.2347267	1.6380612	0.1014089
land_typeSawah	0.4203321	0.1708236	2.4606213	0.0138697
farm_scaleSedang	0.0995129	0.1853786	0.5368091	0.5913995
farm_scaleBesar	1.1680115	0.2896182	4.0329349	0.0000551
farm_scaleSangat Besar	1.3380281	0.4130073	3.2397204	0.0011965
training_count	0.1342151	0.0391168	3.4311406	0.0006010

# Ringkasan model numerik
tidy(model_numeric) %>%
kable(format = "html", booktabs = TRUE, caption = "Ringkasan Koefisien Model dengan Ordinal sebagai Numerik") %>%
kable_styling(latex_options = c("hold_position", "striped"))

## Warning in attr(x, "align"): 'xfun::attr()' is deprecated.
## Use 'xfun::attr2()' instead.
## See help("Deprecated")

Ringkasan Koefisien Model dengan Ordinal sebagai Numerik
term	estimate	std.error	statistic	p.value
(Intercept)	-0.1590135	0.2787187	-0.5705159	0.5683279
land_typeSawah	0.4402752	0.1702842	2.5855310	0.0097229
farm_scale_numeric	0.4582438	0.0970705	4.7207316	0.0000023
training_count	0.1307419	0.0389111	3.3600196	0.0007794

Kesimpulan :

Land_type : Tipe petani Sawah memiliki peluang sukses lebih tinggi dibandingkan petani Ladang.
Farm_Scale :
- Jika dianggap dummy, tiap jenjang dibandingkan dengan skala berukuran Kecil.
- Jika dianggap numeric, tiap kenaikan satu tingkat skala (misal dari kecil ke sedang, sedang ke besar) meningkatkan peluang sukses.
Training_Counts : Jumlah pelatihan yang lebih tinggi meningkatkan peluang sukses.

11 PEMILIHAN MODEL REGRESI LOGISTIK DAN EVALUASI

11.1 MEMBANGUN MODEL REGRESI LOGISTIK : PENDEKATAN CONFIRMATORY DAN EXPLORATORY

Dalam analisis regresi logistik, pemilihan model yang tepat sangat penting untuk memprediksi probabilitas terjadinya suatu peristiwa (dengan respons biner, seperti sukses/gagal). Terdapat dua pendekatan utama dalam membangun model regresi logistik:

Pendekatan Confirmatory (Konfirmatori)

Pendekatan ini digunakan ketika peneliti sudah memiliki hipotesis atau teori yang jelas mengenai hubungan antara variabel prediktor dan variabel respons. Tujuannya bukan sekadar menemukan model terbaik, melainkan menguji kebenaran hubungan yang telah diasumsikan sebelumnya.

Ciri-ciri:
- Model dibangun berdasarkan teori atau penelitian sebelumnya.
- Tujuan utama adalah menguji apakah efek dari variabel tertentu signifikan secara statistik.
- Peneliti sering membangun model penuh terlebih dahulu, lalu menguji apakah penghilangan atau penambahan variabel (misalnya, efek interaksi) meningkatkan kualitas model.
- Uji signifikansi dilakukan dengan membandingkan dua model: model dengan dan tanpa efek tertentu, misalnya melalui Likelihood Ratio Test (LRT).
Contoh: Jika teori menyatakan bahwa variabel x1 dan x2 berpengaruh terhadap kemungkinan seseorang membeli produk, maka model dibangun langsung dengan kedua variabel tersebut. Selanjutnya, diuji apakah kontribusi x2 signifikan dalam menjelaskan probabilitas pembelian.
Pendekatan Exploratory (Eksploratori)

Pendekatan ini digunakan ketika peneliti belum memiliki teori yang pasti atau ingin mengeksplorasi hubungan potensial antara variabel dalam data. Tujuannya adalah menemukan model yang paling efisien dan akurat secara statistik, berdasarkan struktur data yang tersedia.

Ciri-ciri:
- Model dibangun bertahap, dengan proses seleksi berdasarkan kriteria statistik seperti AIC (Akaike Information Criterion), deviance, atau log-likelihood.
- Bertujuan menemukan model yang parsimonious: cukup sederhana, tetapi memiliki performa prediktif yang baik.
- Digunakan teknik seleksi variabel seperti :
  - Forward Selection: dimulai dari model kosong, lalu satu per satu variabel ditambahkan jika signifikan.
  - Backward Elimination: dimulai dari model penuh, lalu variabel yang tidak signifikan secara bertahap dihapus.
  - Stepwise Selection: gabungan dari forward dan backward; variabel bisa masuk dan keluar model secara dinamis.

Tujuan : Menemukan model yang parsimonious, yaitu cukup sederhana namun memiliki performa prediksi yang baik.

Pemilihan antara pendekatan Confrmatory dan Exploratory bergantung pada tujuan penelitian. Jika ingin menguji hipotesis tertentu, gunakan pendekatan Confrmatory. Jika ingin menemukan model terbaik berdasarkan data, gunakan pendekatan Exploratory.

11.2 SIMULASI DATA

Sebuah dataset penjualan terdiri dari beberapa variabel, seperti beli (0 = tidak beli, 1 = beli) sebagai respons, dan umur (x1), gender (x2), dan diskon (x3) sebagai prediktor.

library(knitr)
library(dplyr)
library(ggplot2)
library(MASS)

## 
## Attaching package: 'MASS'

## The following object is masked from 'package:dplyr':
## 
##     select

library(caret)
library(pROC)

## Warning: package 'pROC' was built under R version 4.4.3

## Type 'citation("pROC")' for a citation.

## 
## Attaching package: 'pROC'

## The following objects are masked from 'package:stats':
## 
##     cov, smooth, var

library(DescTools)

## Warning: package 'DescTools' was built under R version 4.4.3

## 
## Attaching package: 'DescTools'

## The following objects are masked from 'package:caret':
## 
##     MAE, RMSE

set.seed(39)
n <- 500
x1 <- rnorm(n)
x2 <- rbinom(n, 1, 0.5)
x3 <- rnorm(n)
lin_pred <- -0.5 + 0.2 * x1 - 0.5 * 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)

##   y         x1 x2         x3
## 1 1 -0.1855657  0 -0.5143321
## 2 0 -1.2292427  0 -0.3466784
## 3 0 -0.4272029  1 -0.5870137
## 4 0 -0.5959819  0 -0.7943017
## 5 0  0.4673236  1  0.4587578
## 6 0  0.4216390  1 -0.6383594

11.3 PEMILIHAN MODEL

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.5695     0.1408  -4.045 5.24e-05 ***
## x1            0.1376     0.1055   1.304    0.192    
## x2           -0.1977     0.2007  -0.985    0.324    
## x3            0.8322     0.1115   7.466 8.27e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 655.68  on 499  degrees of freedom
## Residual deviance: 584.69  on 496  degrees of freedom
## AIC: 592.69
## 
## Number of Fisher Scoring iterations: 3

Interpretasi : Variabel diskon (x3) menunjukkan p-value < 0,05 yang menunjukkan bahwa diskon berpengaruh signifikan dalam meningkatkan peluang beli dari seseorang.

11.4 METODE STEPWISE:FORWARD, BACKWARD, DAN KEDUA ARAH

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

##               df      AIC
## model_full     4 592.6936
## step_forward   2 591.2231
## step_backward  2 591.2231
## step_both      2 591.2231

11.5 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.7206

Interpretasi : Model mampu membedakan pembeli dan non-pembeli sebanyak 72,06%.

11.6 PSEUDO R-SQUARED

PseudoR2(step_both, which = c("CoxSnell", "Nagelkerke", "McFadden"))

##   CoxSnell Nagelkerke   McFadden 
##  0.1279633  0.1751603  0.1044127

Interpretasi : Model yang diperoleh dari ketiga metode pseudo \(R^2-squared\) masih belum terlalu kuat dalam menjelaskan variasi data karena seluruh nilai tersebut masih berada di bawah 0,2 atau model masih menjelaskan kurang dari 20% variasi data. Bisa ditingkatkan dengan menambahkan variabel atau mencari model alternatif.

11.7 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 276 114
##          1  42  68
##                                           
##                Accuracy : 0.688           
##                  95% CI : (0.6454, 0.7284)
##     No Information Rate : 0.636           
##     P-Value [Acc > NIR] : 0.008375        
##                                           
##                   Kappa : 0.2639          
##                                           
##  Mcnemar's Test P-Value : 1.312e-08       
##                                           
##             Sensitivity : 0.3736          
##             Specificity : 0.8679          
##          Pos Pred Value : 0.6182          
##          Neg Pred Value : 0.7077          
##              Prevalence : 0.3640          
##          Detection Rate : 0.1360          
##    Detection Prevalence : 0.2200          
##       Balanced Accuracy : 0.6208          
##                                           
##        'Positive' Class : 1               
##

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

## Sensitivity Specificity 
##   0.3736264   0.8679245

Interpretasi :

Dari seluruh pembeli (sensitivity, aktual 1), hanya 37,36% saja yang berhasil diprediksi oleh model. Sementara itu, dari seluruh non-pembeli (specificity, aktual 0), 86,79% berhasil diprediksi oleh model.
Model lebih baik mengenali non-pembeli dibandingkan pembeli. Hal ini dapat disebabkan oleh adanya ketidakseimbangan data, seperti lebih banyak orang yang tidak membeli dibandingkan membeli.

11.8 METODE PERBANDINGAN MODEL DALAM REGRESI LOGISTIK

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

11.9 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

##     Model      AIC Deviance
## 1 Model 1 656.9178 652.9178
## 2 Model 2 657.2825 651.2825
## 3 Model 3 592.6936 584.6936

Interpretasi : Model 3 memiliki AIC dan deviance paling rendah dibandingkan model 1 dan 2 sehingga model ini merupakan model paling baik di antara ketiganya.

11.10 LIKELIHOOD-RATIO TEST

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

## Analysis of Deviance Table
## 
## Model 1: y ~ x1
## Model 2: y ~ x1 + x2
##   Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1       498     652.92                     
## 2       497     651.28  1   1.6353    0.201

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

## Analysis of Deviance Table
## 
## Model 1: y ~ x1 + x2
## Model 2: y ~ x1 + x2 + x3
##   Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
## 1       497     651.28                          
## 2       496     584.69  1   66.589 3.345e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

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

## Analysis of Deviance Table
## 
## Model 1: y ~ x1
## Model 2: y ~ x1 + x2 + x3
##   Resid. Df Resid. Dev Df Deviance  Pr(>Chi)    
## 1       498     652.92                          
## 2       496     584.69  2   68.224 1.532e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

11.11 PRINSIP PARSIMONY

Model yang kompleks sering memiliki AIC dan deviance yang lebih kecil. Namun, model sederhana lebih mudah diinterpretasikan. Jika penurunan AIC tidak signifikan, pilih model lebih sederhana. Prinsip parsimony mencegah overfitting.

11.11.1 RUMUS DAN PENJELASAN

11.11.1.1 RUMUS AIC

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

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

11.11.1.2 RUMUS DEVIANCE

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

Penjelasan:
Deviance mengukur seberapa jauh model saat ini dibandingkan dengan model sempurna (model saturasi). Nilai deviance yang kecil menunjukkan model memberikan prediksi yang mendekati data aktual.

11.11.1.3 RUMUS LIKELIHOOD-RATIO TEST

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

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

11.12 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, df$y, positive = "1")
conf_matrix

## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   0   1
##          0 272 114
##          1  46  68
##                                           
##                Accuracy : 0.68            
##                  95% CI : (0.6371, 0.7207)
##     No Information Rate : 0.636           
##     P-Value [Acc > NIR] : 0.02209         
##                                           
##                   Kappa : 0.2489          
##                                           
##  Mcnemar's Test P-Value : 1.178e-07       
##                                           
##             Sensitivity : 0.3736          
##             Specificity : 0.8553          
##          Pos Pred Value : 0.5965          
##          Neg Pred Value : 0.7047          
##              Prevalence : 0.3640          
##          Detection Rate : 0.1360          
##    Detection Prevalence : 0.2280          
##       Balanced Accuracy : 0.6145          
##                                           
##        'Positive' Class : 1               
##

Interpretasi :

Model cukup baik dalam mengenali non-pembeli yaitu sebesar 85,5% dan secara keseluruhan mampu menjelaskan sebanyak 68% variasi data. Namun, sensitivitas nya rendah yaitu sebesar 37,36% yang artinya banyak yang pembeli yang terlewatkan dari prediksi model.

11.12.1 SENSITIVITAS DAN SPESIFISITAS

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

\[ \text{Sensitivity} = \frac{TP}{{TP+FN}} \]
Spesifisitas : Kemampuan model mendeteksi kelas negatif secara benar (True Negative Rate)

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

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

## Sensitivity Specificity 
##   0.3736264   0.8553459

11.13 KURVA ROC (RECEIVER OPERATING CHARACTERISTIC)

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

Definisi
- Sumbu Y :
  
  \[\text{Sensitivity} = \frac{TP}{{TP+FN}} \]
- Sumbu X :
  
  \[1-\text{Specificity} = \frac{FP}{{FP+TN}} \]
- Garis diagonal (dari kiri bawah ke kanan atas) menunjukkan performa acak (random guess).
- Kurva yang mendekati pojok kiri atas menunjukkan performa klasifkasi yang lebih baik.
Cut-Off dan Pergerakan Kurva

Saat cut-of menurun, model mengklasifkasikan lebih banyak pengamatan sebagai positif :
- Sensitivitas naik
- Spesifisitas turun
Saat cut-of naik, model menjadi lebih konservatif :
- Sensitivitas turun
- Spesifisitas naik
Kurva ROC Ideal

Kurva ideal memiliki bentuk :
- Naik tajam secara vertikal hingga mencapai sensitivitas = 1
- Lalu bergerak secara horizontal menuju 1 - specifcity = 1
- Area under the curve (AUC) mendekati 1
Interpretasi Luas Area (AUC)
- AUC = 0.5: model tidak lebih baik dari tebak acak
- AUC > 0.7: model cukup baik
- AUC > 0.9: model sangat baik
- AUC dikenal juga sebagai concordance index, yaitu probabilitas bahwa model memberikan nilai skor probabilitas yang lebih tinggi untuk kasus positif daripada kasus negatif.
Kegunaan Kurva ROC
- Untuk membandingkan performa beberapa model klasifkasi
- Untuk memilih threshold (cut-of) optimal berdasarkan kebutuhan aplikasi
Visualisasi dalam R

model <- glm(y ~ x1 + x2 + x3, data = df, family = binomial)
pred <- predict(model, type = "response")
roc_obj <- roc(df$y, pred)

## Setting levels: control = 0, case = 1

## Setting direction: controls < cases

plot(roc_obj)

auc(roc_obj)

## Area under the curve: 0.7224

Simulasi Pemilihan Threshold Optimal

thresholds <- seq(0.1, 0.9, by = 0.05)
results <- data.frame(Threshold = thresholds)

cm <- table(Pred = pred_class, Obs = df$y)
cm

##     Obs
## Pred   0   1
##    0 272 114
##    1  46  68

results$Sensitivity <- sapply(thresholds, function(t) {
pred_class <- ifelse(pred >= t, 1, 0)
cm <- table(factor(pred_class, levels=c(0,1)), factor(df$y, levels=c(0,1)))
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(factor(pred_class, levels=c(0,1)), factor(df$y, levels=c(0,1)))
TN <- cm["0", "0"]
FP <- cm["1", "0"]
TN / (TN + FP)
})
print(results)

##    Threshold Sensitivity Specificity
## 1       0.10 1.000000000  0.03459119
## 2       0.15 0.972527473  0.15094340
## 3       0.20 0.928571429  0.27358491
## 4       0.25 0.868131868  0.42138365
## 5       0.30 0.785714286  0.54716981
## 6       0.35 0.697802198  0.63522013
## 7       0.40 0.593406593  0.74213836
## 8       0.45 0.494505495  0.82075472
## 9       0.50 0.373626374  0.85534591
## 10      0.55 0.318681319  0.90251572
## 11      0.60 0.219780220  0.92138365
## 12      0.65 0.175824176  0.94968553
## 13      0.70 0.065934066  0.97169811
## 14      0.75 0.032967033  0.99371069
## 15      0.80 0.016483516  1.00000000
## 16      0.85 0.005494505  1.00000000
## 17      0.90 0.000000000  1.00000000

Cut-of optimal bisa dipilih berdasarkan maksimum dari Sensitivity + Specificity. Atau mempertimbangkan trade-off sesuai tujuan aplikasi maka prioritaskan sensitivitas tinggi).

Catatan
- ROC cocok saat proporsi kelas seimbang
- Untuk data dengan kelas tidak seimbang, precision-recall curve bisa lebih informatif

11.14 PRECISION-RECALL CURVE (PR CURVE)

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

Definisi
- Precision (Presisi) : Proporsi prediksi positif yang benar-benar positif
  
  \[ \text{Precision}=\frac{TP}{TP+FP} \]
- Recall (Sensitivitas) : Proporsi kasus positif yang berhasil diprediksi positif
  
  \[\text{Recall}=\frac{TP}{TP+FN} \]
Interpretasi
- PR Curve menunjukkan bagaimana presisi berubah saat recall meningkat.
- Idealnya, kita ingin nilai presisi dan recall keduanya tinggi, tetapi biasanya ada trade-of.
- Model dengan performa baik memiliki PR Curve yang melengkung ke pojok kanan atas.
Area Under PR Curve
- Luas kurva (AUPRC) mendekati 1 berarti model sangat baik.
- Baseline AUPRC = prevalensi kelas positif dalam data
PR Curve vs ROC Curve

\[ \begin{array}{|l|l|l|} \hline \textbf{Aspek} & \textbf{ROC Curve} & \textbf{Precision-Recall Curve} \\ \hline \text{Fokus} & \text{Semua kelas} & \text{Kelas positif saja} \\ \text{Kuat di} & \text{Data seimbang} & \text{Data tidak seimbang} \\ \text{Sumbu Y} & \text{Sensitivitas (Recall)} & \text{Precision} \\ \text{Sumbu X} & \text{1 - Spesifisitas} & \text{Recall} \\ \hline \end{array} \]
Visualisasi PR Curve di R

library(PRROC)

## Warning: package 'PRROC' was built under R version 4.4.3

## Loading required package: rlang

## 
## Attaching package: 'rlang'

## The following object is masked from 'package:magrittr':
## 
##     set_names

## The following object is masked from 'package:gtools':
## 
##     chr

model <- glm(y ~ x1 + x2 + x3, data = df, family = binomial)
prob <- predict(model, type = "response")
pr <- pr.curve(scores.class0 = prob[df$y == 1],
scores.class1 = prob[df$y == 0],
curve = TRUE)
plot(pr)

Catatan
- PR Curve sangat informatif untuk aplikasi seperti deteksi penipuan atau diagnosis penyakit langka.
- Gunakan PR Curve saat : – Kelas positif jauh lebih jarang daripada kelas negatif – Tujuan aplikasi lebih mementingkan presisi terhadap kelas minoritas

11.15 PSEUDO R-SQUARED PADA REGRESI LOGISTIK

11.15.1 LIKELIHOOD DAN RUMUS

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

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

Dengan :

\(L_0\) : Likelihood model null

\(L_1\) : Likelihood model penuh

11.15.2 PERHITUNGAN MANUAL R-SQUARED

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

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

##   R2_Cox_Snell R2_McFadden
## 1    0.1323638   0.1082706

11.15.3 PERHITUNGAN OTOMATIS R-SQUARED

# Menggunakan pscl
if (!require(pscl)) install.packages("pscl"); library(pscl)

## Loading required package: pscl

## Warning: package 'pscl' was built under R version 4.4.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 
## -292.3467850 -327.8424924   70.9914149    0.1082706    0.1323638    0.1811840

# Menggunakan rcompanion
if (!require(rcompanion)) install.packages("rcompanion"); library(rcompanion)

## Loading required package: rcompanion

## Warning: package 'rcompanion' was built under R version 4.4.3

nagelkerke(model)

## $Models
##                                             
## Model: "glm, y ~ x1 + x2 + x3, binomial, df"
## Null:  "glm, y ~ 1, binomial, df"           
## 
## $Pseudo.R.squared.for.model.vs.null
##                              Pseudo.R.squared
## McFadden                             0.108271
## Cox and Snell (ML)                   0.132364
## Nagelkerke (Cragg and Uhler)         0.181184
## 
## $Likelihood.ratio.test
##  Df.diff LogLik.diff  Chisq    p.value
##       -3     -35.496 70.991 2.6179e-15
## 
## $Number.of.observations
##           
## Model: 500
## Null:  500
## 
## $Messages
## [1] "Note: For models fit with REML, these statistics are based on refitting with ML"
## 
## $Warnings
## [1] "None"

# Menggunakan DecTools
if (!require(DescTools)) install.packages("DescTools"); library(DescTools)
PseudoR2(model, which = "all")

##        McFadden     McFaddenAdj        CoxSnell      Nagelkerke   AldrichNelson 
##      0.10827061      0.09606963      0.13236385      0.18118399      0.12433009 
## VeallZimmermann           Efron McKelveyZavoina            Tjur             AIC 
##      0.21913941      0.13363464      0.18396688      0.13507311    592.69356996 
##             BIC          logLik         logLik0              G2 
##    609.55200235   -292.34678498   -327.84249244     70.99141492

12 DISTRIBUSI MULTINOMIAL

Distribusi multinomial adalah perluasan dari distribusi binomial untuk lebih dari dua kategori. Jika X1, X2, …, X𝑘 menyatakan banyaknya kejadian dalam masing-masing dari 𝑘 kategori, maka:

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

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

12.1 STUDI KASUS

Sebuah perusahaan melakukan survei terhadap 20 karyawan tentang moda transportasi yang digunakan ke kantor: Motor (M), Mobil (Mo), Kereta (K)

Hasil survei:

Motor : 10 orang
Mobil : 6 orang
Kereta : 4 orang

Probabilitas teoretik yang diharapkan:

\(p_{M}\)=0,4
\(p_{Mo}\)=0,3
\(p_{K}\)=0,3

Pertanyaan : Berapa peluang bahwa hasil survei menunjukkan 10 orang naik motor, 6 orang naik mobil, dan 4 orang naik kereta?

Rumus distribusi multinomial :

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

dimana :

\(n\) = 20, \(x_{1}\) = 10, \(x_{2}\) = 6, \(x_{3}\) = 4
\(p_{1}\) = 0,4, \(p_{2}\) = 0,3, \(p_{3}\) = 0,3

Perhitungan Manual di R

n <- 20
x <- c(10, 6, 4)
p <- c(0.4, 0.3, 0.3)

# Hitung komponen-koefisien
faktorial_total <- factorial(n)
faktorial_x <- prod(factorial(x))
koefisien <- faktorial_total / faktorial_x
# Hitung peluang
peluang <- koefisien * prod(p^x)
peluang

## [1] 0.02402317

Interpretasi : Peluang bahwa hasil survei menunjukkan 10 orang naik motor, 6 orang naik mobil, dan 4 orang naik kereta adalah sebesar 0,025.

12.2 REGRESI LOGISTIK MULTINOMIAL

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

Misalkan Y memiliki K kategori, dan kita pilih baseline kategori K, maka model logit untuk kategori j adalah sebagai berikut. \[ \log \left( \frac{P(Y = j)}{P(Y = K)} \right) = \beta_{j0} + \beta_{j1} x_1 + \cdots + \beta_{jp} x_p \] untuk \(( j = 1, 2, \ldots, K-1 )\).

12.2.1 BASELINE-CATEGORY LOGIT MODEL

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

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

dengan: - \((\pi_j)\) adalah probabilitas respon berada di kategori j - \((\pi_c)\) adalah probabilitas respon berada di kategori baseline Maka, terdapat sebanyak \((c-1)\) fungsi logit.

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

12.3 MODEL REGRESI

Jika terdapat satu prediktor x, maka bentuk umum model logit-nya adalah: \[ \log \left( \frac{\pi_j}{\pi_c} \right) = \alpha_j + \beta_j x, \quad j = 1, \ldots, c-1 \]

12.4 CONTOH KASUS : 3 KATEGORI RESPON

Misalkan respon Y memiliki tiga kategori: \(Y \in \{1, 2, 3\}\), dan kita gunakan kategori ke-3 sebagai baseline. Maka:

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

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

12.5 RELASI ANTAR KATEGORI

Jika ingin menghitung logit antara kategori 1 dan 2: \[ \begin{aligned} \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 \end{aligned} \]

12.6 KARAKTERISTIK MODEL BASELINE-CATEGORY LOGIT

Digunakan untuk respon dengan kategori >2
Menghasilkan \((c-1)\) fungsi logit terhadap satu baseline
Logit antar kategori selain baseline dapat dihitung dari selisih dua logit terhadap baseline
Implementasi di R menggunakan fungsi multinom() dari package nnet, dan kategori baseline bisa diatur dengan relevel()

12.7 ESTIMASI PARAMETER

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

Fungsi log-likelihood: \[ \ell(\beta) = \sum_{i=1}^{n} \sum_{j=1}^{K} y_{ij} \log(\pi_{ij}) \] dengan: - \(\pi_{ij} = P(Y_i = j | x_i)\) - \(y_{ij} = 1\) jika \(Y_i = j\), 0 jika tidak.

12.8 CONTOH KASUS

Sebuah perusahaan pengembang aplikasi ingin memahami faktor-faktor yang memengaruhi frekuensi penggunaan aplikasi oleh pengguna. Frekuensi dikategorikan sebagai : Rendah, Sedang, dan Tinggi.

Perusahaan mengumpulkan data dari 150 pengguna aplikasi dengan informasi berikut:

- Platform: Jenis sistem operasi yang digunakan (Android, iOS)

- Age: Usia pengguna

- Occupation: Pekerjaan pengguna (Pelajar, Karyawan, Wirausaha)

- DailyTime: Rata-rata waktu penggunaan aplikasi per hari (dalam jam)

12.9 SIMULASI DATA

set.seed(123)
n <- 150

Platform <- sample(c("Android", "iOS"), n, replace = TRUE)
Age <- round(rnorm(n, mean = 35, sd = 8))
Occupation <- sample(c("Pelajar", "Karyawan", "Wirausaha"), n, replace = TRUE)
DailyTime <- round(rnorm(n, mean = 9, sd = 3), 0)

# Simulasi Frekuensi berdasarkan Occupation
Frequency <- sapply(Occupation, function(job) {
  if (job == "Pelajar") {
    sample(c("Rendah", "Sedang", "Tinggi"), size = 1, prob = c(0.6, 0.2, 0.2))
  } else if (job == "Karyawan") {
    sample(c("Rendah", "Sedang", "Tinggi"), size = 1, prob = c(0.3, 0.3, 0.4))
  } else {
    sample(c("Rendah", "Sedang", "Tinggi"), size = 1, prob = c(0.4, 0.4, 0.2))
  }
})

df <- data.frame(
  Frequency = factor(Frequency),
  Platform = factor(Platform),
  Age = Age,
  Occupation = factor(Occupation),
  DailyTime = DailyTime
)

# Set baseline kategori
df$Frequency <- relevel(df$Frequency, ref = "Sedang")
head(df)

##   Frequency Platform Age Occupation DailyTime
## 1    Tinggi  Android  43   Karyawan        11
## 2    Sedang  Android  33    Pelajar         5
## 3    Sedang  Android  25   Karyawan         6
## 4    Rendah      iOS  36    Pelajar         6
## 5    Tinggi  Android  34   Karyawan         9
## 6    Tinggi      iOS  35    Pelajar         7

13 REGRESI LOGISTIK ORDINAL

Regresi logistik ordinal digunakan ketika variabel respon Y bersifat ordinal (memiliki urutan), misalnya tingkat kepuasan : Rendah, Sedang, Tinggi.

Model ini berbeda dengan :

Regresi logistik biner : hanya 2 kategori
Regresi logistik multinomial : kategori > 2 tetapi tidak berurutan

13.1 KONSEP CUMULATIVE LOGIT MODEL

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

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

\(\alpha_{j}\) : intersep khusus untuk kategori ke-j
\(\beta\) : koefisien regresi (sama untuk semua kategori kumulatif)

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

13.2 INTERPRETASI KOEFISIEN

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

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

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

Odds ratio : \[OR=e^{\beta}\]

13.3 SIMULASI DATA

Sebuah studi kesehatan ingin mengetahui apakah jumlah aktivitas yang dilakukan dalam satu minggu memengaruhi tingkat stres individu.

Tingkat stres diklasifikasikan menjadi : Rendah, Sedang, Tinggi.

Diperoleh data dari 160 responden dengan variabel:

Aktivitas : Jumlah aktivitas per minggu
Stres : Tingkat stres (ordinal)

set.seed(100)
n <- 160 
aktivitas <- round(runif(n, 0, 12))  
# jam per minggu 
stres <- cut(   2.5 + 0.5*aktivitas + rnorm(n),   breaks = c(-Inf, 1.5, 3, Inf),   labels = c("Rendah", "Sedang", "Tinggi"),   ordered_result = TRUE ) 
df <- data.frame(stres, aktivitas) 
head(df)

##    stres aktivitas
## 1 Tinggi         4
## 2 Tinggi         3
## 3 Tinggi         7
## 4 Rendah         1
## 5 Tinggi         6
## 6 Tinggi         6

table(df$stres)

## 
## Rendah Sedang Tinggi 
##      3     14    143

13.4 ESTIMASI MODEL ORDINAL

library(MASS) 
model_ord <- polr(stres ~ aktivitas, data = df, Hess = TRUE) 
summary(model_ord)

## Call:
## polr(formula = stres ~ aktivitas, data = df, Hess = TRUE)
## 
## Coefficients:
##            Value Std. Error t value
## aktivitas 0.7756     0.1806   4.295
## 
## Intercepts:
##               Value   Std. Error t value
## Rendah|Sedang -1.6785  0.6616    -2.5371
## Sedang|Tinggi  0.5860  0.4868     1.2037
## 
## Residual Deviance: 82.24959 
## AIC: 88.24959

13.5 NILAI P-VALUE

(ctable <- coef(summary(model_ord)))

##                    Value Std. Error   t value
## aktivitas      0.7755836  0.1805664  4.295282
## Rendah|Sedang -1.6785056  0.6615848 -2.537098
## Sedang|Tinggi  0.5859999  0.4868323  1.203700

p <- pnorm(abs(ctable[, "t value"]), lower.tail = FALSE) * 2 
ctable <- cbind(ctable, "p value" = round(p, 2))
ctable

##                    Value Std. Error   t value p value
## aktivitas      0.7755836  0.1805664  4.295282    0.00
## Rendah|Sedang -1.6785056  0.6615848 -2.537098    0.01
## Sedang|Tinggi  0.5859999  0.4868323  1.203700    0.23

Interpretasi :

Jumlah aktivitas berpengaruh signifikan terhadap tingkat stress karena p-value < 0,05. Koefisien yang positif juga menandakan bahwa semakin banyak aktivitas yang diikuti selama satu minggu maka semakin besar kemungkinan seseorang untuk mengalami tingkat stres yang lebih tinggi.
Cut-point menanadakan batas logit antara kategori stres. P-value yang signifikan hanya untuk cut-point pertama menunjukkan bahwa model lebih dapat membedakan antara “Rendah” dan “Sedang”, tetapi tidak untuk “Sedang” dan”Tinggi”.

13.6 PREDIKSI PROBABILITAS

newdata <- data.frame(aktivitas = 5:9) 
predict(model_ord, newdata = newdata, type = "probs")

##         Rendah      Sedang    Tinggi
## 1 0.0038477043 0.032001788 0.9641505
## 2 0.0017753034 0.015056640 0.9831681
## 3 0.0008181956 0.007002851 0.9921790
## 4 0.0003768923 0.003239452 0.9963837
## 5 0.0001735697 0.001494777 0.9983317

13.7 GOODNESS-OF-FIT DAN PROPORTIONAL ODDS

Model cumulative logit mengasumsikan bahwa pengaruh setiap prediktor bersifat konstan di seluruh batas kategori respons (asumsi proportional odds).

13.8 ALTERNATIF MODEL ORDINAL

Selain cumulative logit, model ordinal lainnya adalah :

Adjacent-category logit
Continuation-ratio (sequential) logit

Model tersebut dapat digunakan saat asumsi proportional odds tidak terpenuhi.

13.9 ASUMSI PARALELISME DALAM REGRESI LOGISTIK ORDINAL

Model 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).

13.9.1 DEFINISI ASUMSI PARALELISME

Asumsi paralelisme menyatakan bahwa koefisien regresi sama untuk setiap kategori kumulatif dari variabel respon.

Bentuk umum model :

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

untuk \(j=1, ..., c-1\)

Hanya intersep (\[\alpha_{j}\]) yang berbeda-beda.
Koefisien \[\beta\] tetap sama untuk semua fungsi logit kumulatif.

13.9.2 VISUALISASI

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

13.9.3 KONSEKUENSI PELANGGARAN ASUMSI

Jika asumsi ini tidak terpenuhi :

Efek prediktor berbeda untuk setiap batas kategori
Model cumulative logit tidak valid
Perlu gunakan model alternatif :
- Generalized Ordinal Logistic Regression
- Partial Proportional Odds Model

13.9.4 PENGUJIAN ASUMSI PARALELISME

Untuk memeriksa validitas asumsi, dapat digunakan :

Likelihood Ratio Test antara model proportional dan non-proportional
Brant test
Jika p-value < 0,05 maka asumsi tidak terpenuhi

library(brant)

## Warning: package 'brant' was built under R version 4.4.3

brant(model_ord)

## -------------------------------------------- 
## Test for X2  df  probability 
## -------------------------------------------- 
## Omnibus      0.97    1   0.32
## aktivitas    0.97    1   0.32
## -------------------------------------------- 
## 
## H0: Parallel Regression Assumption holds

14 LOG LINEAR MODEL

Analisis data kategorik merupakan bagian penting dalam statistika terapan karena banyak fenomena di dunia nyata yang menghasilkan data dalam bentuk kategori, seperti jenis kelamin, status pekerjaan, tingkat pendidikan, preferensi konsumen, atau hasil diagnosis medis. Data kategori ini umumnya dianalisis menggunakan tabel kontingensi, model log-linier, dan model regresi logistik. Masing-masing pendekatan memiliki kekuatan dan kelemahan tergantung pada tujuan analisis dan struktur data.

Tabel kontingensi digunakan sebagai langkah awal eksplorasi untuk melihat hubungan antara dua atau lebih variabel kategorik. Misalnya, dalam studi tentang efek obat terhadap serangan jantung, tabel kontingensi dapat menyajikan jumlah pasien yang mengalami atau tidak mengalami serangan jantung, berdasarkan jenis obat yang dikonsumsi. Tabel ini membantu mengidentifikasi pola awal dan menghitung ukuran asosiasi seperti odds ratio, risk ratio, dan chi-square statistic untuk menguji independensi antar variabel.

Namun, ketika ingin membangun model statistik yang dapat mengendalikan efek dari banyak variabel dan interaksinya secara simultan, maka model log-linier menjadi sangat berguna. Model log-linier adalah bentuk khusus dari Generalized Linear Model (GLM) yang digunakan pada frekuensi sel dalam tabel kontingensi dan mengasumsikan distribusi Poisson. Tidak seperti regresi logistik, model log-linier tidak menetapkan variabel mana yang dependen dan mana yang independen, karena seluruh variabel diperlakukan secara simetris. Model ini lebih cocok ketika tujuan analisis adalah untuk memahami struktur asosiasi atau independensi antar variabel, bukan untuk prediksi.

Struktur model log-linier ditentukan berdasarkan efek utama dari masing-masing variabel serta interaksi di antara variabel-variabel tersebut. Misalnya, dalam tabel kontingensi tiga arah (misalnya: jenis kelamin, status merokok, dan penyakit paru), model log-linier dapat membedakan apakah interaksi dua variabel cukup menjelaskan data, ataukah diperlukan interaksi tiga arah untuk menangkap struktur asosiasinya. Penyesuaian model dapat dilakukan menggunakan metode likelihood ratio test untuk membandingkan model sederhana dengan model lebih kompleks.

Di sisi lain, regresi logistik adalah pendekatan paling umum ketika terdapat satu variabel kategorik yang secara eksplisit dianggap sebagai variabel dependen (misalnya, kejadian penyakit: ya/tidak), dan satu atau lebih variabel kategorik atau numerik sebagai prediktor. Model ini memodelkan logit dari probabilitas kejadian (yaitu log odds), dan sangat berguna dalam studi observasional dan eksperimental untuk menjelaskan atau memprediksi peluang suatu outcome. Regresi logistik juga memiliki ekstensi untuk outcome kategorik lebih dari dua kelas, seperti regresi logistik multinomial dan regresi logistik ordinal.

Dengan demikian, meskipun ketiganya beroperasi pada data kategorik, tabel kontingensi bersifat deskriptif, model log-linier bersifat eksploratif terhadap hubungan simetris, sedangkan regresi logistik bersifat prediktif terhadap outcome kategorik. Pemilihan metode tergantung pada apakah fokus utama analisis adalah deskripsi, eksplorasi struktur, atau prediksi hasil berdasarkan variabel penjelas. Kombinasi dari ketiganya sering digunakan dalam praktik untuk memperoleh pemahaman komprehensif dari data kategorik yang dianalisis.

Ringkasan Dalam analisis data kategorik, terdapat beberapa pendekatan statistik yang umum digunakan, antara lain:

Tabel Kontingensi: penyajian frekuensi gabungan dari dua atau lebih variabel kategorik.
Model Loglinear: digunakan untuk memodelkan struktur asosiasi di dalam tabel kontingensi tanpa menganggap ada variabel dependen.
Model Regresi Logistik: digunakan untuk memodelkan probabilitas dari kategori variabel dependen berdasarkan variabel independen.

Meskipun ketiganya dapat digunakan pada data kategorik, pendekatan dan interpretasinya sangat berbeda.

Tabel Kontingensi

Tabel kontingensi menyajikan jumlah frekuensi dari kombinasi kategori antar variabel.

Contoh:

Sebuah studi dilakukan untuk menguji apakah program relaksasi mingguan dapat mengurangi kejadian serangan cemas pada penderita gangguan kecemasan umum (GAD). Sebanyak 100 partisipan dibagi ke dalam dua kelompok:

Kelompok intervensi mengikuti program relaksasi terstruktur.
Kelompok kontrol tidak mengikuti program, hanya diberi edukasi umum (placebo).

Setelah 4 minggu, peneliti mencatat jumlah peserta di setiap kelompok yang masih mengalami serangan cemas dan yang tidak dalam tabel berikut.

table_data <- matrix(c(55, 35, 70, 30), nrow=2, 
       dimnames = list(Obat = c("Program_Relaksasi", "Placebo"),
                       Serangan = c("Ya", "Tidak")))
table_data

##                    Serangan
## Obat                Ya Tidak
##   Program_Relaksasi 55    70
##   Placebo           35    30

Tabel kontingensi bersifat deskriptif dan tidak melibatkan pemodelan probabilitas.

Model Loglinear

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

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

Cocok untuk menyelidiki asosiasi dan independensi antar variabel.
Tidak membedakan antara variabel respon dan penjelas.
Umumnya digunakan pada tabel multi-dimensi (2 arah, 3 arah, dst).

library(MASS)
loglm(~Obat*Serangan, data = table_data)

## Call:
## loglm(formula = ~Obat * Serangan, data = table_data)
## 
## Statistics:
##                  X^2 df P(> X^2)
## Likelihood Ratio   0  0        1
## Pearson            0  0        1

Interpretasi :

Tidak ada eror atau tidak ada penyimpangan antara data yang diamati dan yang diprediksi oleh model.
Derajat bebas (df) = 0, artinya tidak ada sisa variasi yang bisa dianalisis.

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:

Sebuah penelitian dilakukan untuk mengevaluasi efektivitas dua jenis program diet: Diet Keto dan Diet Mediterania terhadap keberhasilan penurunan berat badan dalam 3 bulan. Keberhasilan didefinisikan sebagai penurunan berat badan ≥ 5 kg. Hasil penelitian terhadap 120 orang adalah sebagai berikut :

# Data frame penelitian
data_diet <- data.frame(
  Sukses = c(1, 0, 1, 0),
  Diet = factor(c("Keto", "Keto", "Mediterania", "Mediterania")),
  Frek = c(40, 20, 30, 30)
)

# Model regresi logistik
model_diet <- glm(Sukses ~ Diet, weights = Frek, family = binomial, data = data_diet)

# Ringkasan model
summary(model_diet)

## 
## Call:
## glm(formula = Sukses ~ Diet, family = binomial, data = data_diet, 
##     weights = Frek)
## 
## Coefficients:
##                 Estimate Std. Error z value Pr(>|z|)  
## (Intercept)       0.6931     0.2739   2.531   0.0114 *
## DietMediterania  -0.6931     0.3764  -1.842   0.0655 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 163.01  on 3  degrees of freedom
## Residual deviance: 159.56  on 2  degrees of freedom
## AIC: 163.56
## 
## Number of Fisher Scoring iterations: 4

Interpretasi :

Intercept menunjukkan log-odds keberhasilan untuk Diet Keto (baseline). Karena p-value nya bernilai 0.0114 < 0.05, maka Diet Keto berhasil menurunkan berat badan dalam 3 bulan.
Koefisien DietMediterania = -0.6931 berarti Diet Mediterania menurunkan log-odds keberhasilan dibanding Diet Keto.
p-value = 0.0655 > 0.05 → berarti perbedaan antara kedua diet tidak signifikan secara statistik. Artinya, jenis diet tidak berpengaruh secara signifikan terhadap keberhasilan penurunan berat badan.

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

14.1 TABEL KONTINGENSI DAN MODEL LOGLINEAR

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

# Contoh tabel 2x2
matrix(c(30, 120, 85, 35), nrow=2,
       dimnames = list(Obat = c("Timolol", "Placebo"),
                       Serangan = c("Ya", "Tidak")))

##          Serangan
## Obat       Ya Tidak
##   Timolol  30    85
##   Placebo 120    35

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

14.2 MODEL SATURATED

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

Cocok sempurna terhadap data
Tidak mengasumsikan independensi antar variabel

Contoh formulasi untuk tabel 2x2:

# Data
library(MASS)
data <- matrix(c(30, 85, 120, 35), nrow=2, byrow=TRUE)
dimnames(data) <- list(Obat = c("Timolol", "Placebo"), Serangan = c("Ya", "Tidak"))
ftable(data)

##         Serangan  Ya Tidak
## Obat                      
## Timolol           30    85
## Placebo          120    35

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

model_saturated <- loglm(~ Obat * Serangan, data = data)
summary(model_saturated)

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

Interpretasi :

Model ini menggunakan semua informasi dalam tabel kontingensi (termasuk interaksi). Akibatnya, tidak ada sisa variasi untuk diuji karena model memperkirakan setiap sel persis seperti data aslinya.
Oleh karena itu, deviasi = 0, dan derajat bebas (df) = 0. Artinya, model prediksi sama hasilnya dengan data asli sehingga model ini fit secara sempurna.

14.3 MODEL INDEPENDENT

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

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

model_indep <- loglm(~ Obat + Serangan, data = data)
summary(model_indep)

## Formula:
## ~Obat + Serangan
## attr(,"variables")
## list(Obat, Serangan)
## attr(,"factors")
##          Obat Serangan
## Obat        1        0
## Serangan    0        1
## attr(,"term.labels")
## [1] "Obat"     "Serangan"
## attr(,"order")
## [1] 1 1
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
## 
## Statistics:
##                       X^2 df P(> X^2)
## Likelihood Ratio 73.35801  1        0
## Pearson          70.45372  1        0

Interpretasi :

Terdapat deviasi besar antara model dan data, yaitu sebesar 73.36.
Model ini adalah model tanpa interaksi, artinya model mengasumsikan bahwa obat dan serangan saling independen. Dan karena p-value < 0.05, maka H0 ditolak. Artinya, terdapat hubungan antara obat dan serangan.
Model independen tidak cocok dengan data, sehingga ada hubungan signifikan antara jenis obat dan jumlah serangan. Artinya, jenis obat mempengaruhi serangan.

14.4 ODDS RATIO

Odds ratio untuk tabel 2x2:

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

Interpretasi nilai OR:

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

14.5 ESTIMASI PARAMETER

Dalam model saturated:

Estimasi dilakukan dengan pembatasan seperti sum-to-zero
Estimasi parameter dilakukan dengan iterative proportional fitting (IPF)

# Estimasi odds ratio dan log-odds
logOR <- log((data[1,1] * data[2,2]) / (data[1,2] * data[2,1]))
logOR

## [1] -2.273598

OR <- exp(logOR)
OR

## [1] 0.1029412

Interpretasi : Odds serangan di grup Placebo adalah 10.3x lebih besar dibanding grup Timolol.

14.6 MODEL LEBIH SEDERHANA DAN PERBANDINGAN MODEL

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

anova(model_indep, model_saturated)

## LR tests for hierarchical log-linear models
## 
## Model 1:
##  ~Obat + Serangan 
## Model 2:
##  ~Obat * Serangan 
## 
##           Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1   73.35801  1                                    
## Model 2    0.00000  0   73.35801         1              0
## Saturated  0.00000  0    0.00000         0              1

Interpretasi : Berdasarkan uji Likelihood Ratio, terdapat perbedaan signifikan antara model independen dan interaksi (delta(dev) = 73.36, df=1, dan p-value < 0.05). Hal ini menunjukkan bahwa terdapat hubungan yang signifikan antara jenis obat dan kejadian serangan.

14.7 STUDI KASUS

Sebuah lembaga riset pendidikan ingin mengetahui apakah terdapat hubungan antara tingkat pendidikan orang tua dan kemungkinan anak mengikuti les privat di luar sekolah. Survei dilakukan terhadap 800 siswa SMA, dan hasilnya diklasifikasikan sebagai berikut :

data_pendidikan <- matrix(c(65, 85,
                            120, 180,
                            90, 160),
                          nrow = 3, byrow = TRUE,
                          dimnames = list(Pendidikan_Ortu = c("SD", "SMA", "Sarjana"),
                                          Les = c("Tidak", "Ya")))
ftable(data_pendidikan)

##                 Les Tidak  Ya
## Pendidikan_Ortu              
## SD                     65  85
## SMA                   120 180
## Sarjana                90 160

loglm(~ Pendidikan_Ortu + Les, data = data_pendidikan)

## Call:
## loglm(formula = ~Pendidikan_Ortu + Les, data = data_pendidikan)
## 
## Statistics:
##                       X^2 df  P(> X^2)
## Likelihood Ratio 2.226945  2 0.3284166
## Pearson          2.226025  2 0.3285677

Interpretasi :

Model yang diujikan adalah model tanpa interaksi. Artinya, model ini mengasumsikan independensi antara pendidikan orang tua dan keputusan anak ikut les.
Berdasarkan hasil Likelihood Ratio (\(G^2\) = 2.23, df = 2, dan p > 0.05), maka H0 diterima. Artinya, pendidikan orang tua dan keputusan anak ikut les tidak berhubungan atau saling independen.

14.8 RANGKUMAN MODEL LOG LINEAR 2 ARAH

14.9 MODEL LOGLINEAR PADA TABEL KONTINGENSI

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

Model log-linear digunakan untuk memodelkan frekuensi (count) pada tabel kontingensi dan menguji asosiasi antar variabel kategorik, tanpa menganggap ada variabel respon dan prediktor.
Model regresi logistik digunakan untuk memodelkan probabilitas kejadian suatu outcome (biner) berdasarkan satu atau lebih prediktor (bisa kategorik maupun kontinu).

14.10 ESTIMASI PARAMETER LOG LINEAR 2 ARAH

Sistem Persamaan Model Log-Linear

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

Constraint Sum-to-Zero

\[ \lambda^A_1 + \lambda^A_2 = 0 \]

\[ \lambda^B_1 + \lambda^B_2 = 0 \]

\[ \lambda^{AB}_{11} + \lambda^{AB}_{12} + \lambda^{AB}_{21} + \lambda^{AB}_{22} = 0 \]

Rumus Estimasi Parameter dengan Sum-to-Zero Constraint

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

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

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

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

Sebuah studi kesehatan dilakukan untuk meneliti apakah ada hubungan antara kebiasaan sarapan pagi dan kejadian kelelahan saat belajar di sekolah. Sebanyak 100 siswa SMA diobservasi dan diklasifikasikan berdasarkan apakah mereka sarapan setiap pagi dan apakah mereka sering merasa lelah saat belajar di kelas. Berikut hasil pengamatan :

Sarapan	Lelah	Tidak Lelah
Ya	25	35
Tidak	30	10

14.11 BENTUK MODEL LOG-LINEAR

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

14.12 ESTIMASI PARAMETER MODEL (MANUAL, SUM-TO-ZERO)

Misalkan:

A1 = Sarapan (Ya), A2 = Tidak
B1 = Lelah, B2 = Tidak lelah

Observasi:

\[ \begin{aligned} n_{11} &= 25,\quad n_{12} = 35 \\ n_{21} &= 30,\quad n_{22} = 10 \end{aligned} \]

\[ \log(\mu_{11}) = \lambda + \lambda^A_1 + \lambda^B_1 + \lambda^{AB}_{11} \]

\[ \log(\mu_{12}) = \lambda + \lambda^A_1 + \lambda^B_2 + \lambda^{AB}_{12} \]

\[ \log(\mu_{21}) = \lambda + \lambda^A_2 + \lambda^B_1 + \lambda^{AB}_{21} \]

\[ \log(\mu_{22}) = \lambda + \lambda^A_2 + \lambda^B_2 + \lambda^{AB}_{22} \]

Constraint Sum-to-Zero

\[ \lambda^A_1 + \lambda^A_2 = 0 \]

\[ \lambda^B_1 + \lambda^B_2 = 0 \]

Langkah-langkah:

Hitung rata-rata log frekuensi sel:

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

\[ = \frac{1}{4} \left[ \log(25) + \log(35) + \log(30) + \log(10) \right] \]

\[ = 3.1195 \]

Efek utama A (Sarapan):

\[\lambda_1^A = \frac{1}{2} \left[ (\log(25) + \log(35)) - (\log(30) + \log(10)) \right] \]

\[ = 0.5352 \]

\[ \lambda_2^A = -0.5352 \]

Efek utama B (Kelelahan):

\[ \lambda_{1}^{B} = \frac{1}{2}\left[(\log(25)+\log(30))-(\log(35)+\log(10))\right] \]

\[ = 0.3811 \]

\[ \lambda_{2}^{B} = -0.3811 \]

Efek interaksi:

\[ \lambda_{11}^{AB} = \frac{1}{4}\left[\log(25)-\log(35)-\log(30)+\log(10)\right] \]

\[ = 0.3587 \]

Ringkasan Parameter

\[ \lambda_{11}^{AB} = \lambda_{22}^{AB}=-0.3587 \]

\[ \lambda_{12}^{AB} = \lambda_{21}^{AB}=+0.3587 \]

14.13 HITUNG ODDS RATIO DAN INTERVAL

\[ OR = \frac{n_{11} n_{22}}{n_{12} n_{21}} = \frac{25 \times 10}{35 \times 30} = \frac{250}{1050} = 0.24 \]

Log odds ratio: \[ \log(OR) = \log(0.24) = -1.4272 \]

Standard error dihitung sebagai:

\(SE = \sqrt{ \frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}} }\)

Dengan substitusi:

\(SE = \sqrt{ \frac{1}{25} + \frac{1}{35} + \frac{1}{30} + \frac{1}{10} } = \sqrt{0.201} = 0.4484\) 95% Confidence Interval for \(\log(OR)\):

\[ \log(OR) \pm 1.96 \times SE = -1.4272 \pm 1.96 \times 0.4484 \]

\[ = (-1.4272 - 0.879), \, (-1.4272 + 0.879) \]

\[ = (-2.3062, -0.5482) \]

Back-transform to get CI for \(OR\):

Lower = \(\exp(-2.3062) = 0.0986\)

Upper = \(\exp(-0.5482) = 0.5741\)

Jadi, OR=0.24 dengan taraf signifikansi 95% memiliki interval kepercayaan dalam rentang 0.0986 - 0.5741.

14.14 FITTING MODEL LOG-LINEAR DENGAN R

# Data 2x2
tabel <- matrix(c(25, 35, 30, 10), nrow = 2, byrow = TRUE)
colnames(tabel) <- c("Lelah", "Tidak lelah")
rownames(tabel) <- c("Sarapan", "Tidak sarapan")
tabel

##               Lelah Tidak lelah
## Sarapan          25          35
## Tidak sarapan    30          10

#Ubah menjadi bentuk tabel untuk glm
data <- as.data.frame(as.table(tabel))
colnames(data) <- c("Sarapan", "Status", "Freq")
data

##         Sarapan      Status Freq
## 1       Sarapan       Lelah   25
## 2 Tidak sarapan       Lelah   30
## 3       Sarapan Tidak lelah   35
## 4 Tidak sarapan Tidak lelah   10

# Model tanpa interaksi
fit_no_inter <- glm(Freq ~ Sarapan + Status, family = poisson, data = data)
summary(fit_no_inter)

## 
## Call:
## glm(formula = Freq ~ Sarapan + Status, family = poisson, data = data)
## 
## Coefficients:
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)            3.4965     0.1576  22.181   <2e-16 ***
## SarapanTidak sarapan  -0.4055     0.2041  -1.986    0.047 *  
## StatusTidak lelah     -0.2007     0.2010  -0.998    0.318    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 16.167  on 3  degrees of freedom
## Residual deviance: 11.138  on 1  degrees of freedom
## AIC: 37.001
## 
## Number of Fisher Scoring iterations: 4

# Model dengan interaksi
fit_inter <- glm(Freq ~ Sarapan * Status, family = poisson, data = data)
summary(fit_inter)

## 
## Call:
## glm(formula = Freq ~ Sarapan * Status, family = poisson, data = data)
## 
## Coefficients:
##                                        Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                              3.2189     0.2000  16.094   <2e-16 ***
## SarapanTidak sarapan                     0.1823     0.2708   0.673   0.5008    
## StatusTidak lelah                        0.3365     0.2619   1.285   0.1988    
## SarapanTidak sarapan:StatusTidak lelah  -1.4351     0.4493  -3.194   0.0014 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 1.6167e+01  on 3  degrees of freedom
## Residual deviance: 2.4425e-15  on 0  degrees of freedom
## AIC: 27.863
## 
## Number of Fisher Scoring iterations: 3

14.15 INTERPRETASI PARAMETER

Berdasarkan model independensi, terdapat pengaruh dari sarapan terhadap status kelelahan. Namun, tidak ada efek signifikan dari status kelelahan secara independen.
Berdasarkan model saturated, pengaruh utama (sarapan dan status kelelahan) tidak berpengaruh secara signifikan, tetapi pengaruhnya sangat signifikan. Artinya, pengaruh sarapan terhadap status kelelahan bergantung pada interaksi antara keduanya.

14.16 ANALISIS DATA TABEL KONTINGENSI 2X3

Sebuah penelitian dilakukan untuk mengetahui hubungan antara Jenis Pekerjaan (Kantoran/Lapangan) dan Tingkat Stres (Rendah/Sedang/Tinggi) pada sekelompok karyawan.

	Rendah	Sedang	Tinggi
Kantoran	14	18	10
Lapangan	12	20	18

14.17 BENTUK MODEL LOG-LINEAR UNTUK TABEL 2X3

Bentuk umum model log-linear untuk tabel 2x3 (dengan sum-to-zero constraint): \[ \log(\mu_{ij}) = \lambda + \lambda_i^{A} + \lambda_j^{B} + \lambda_{ij}^{AB} \] dengan:

\(\mu_{ij}\): ekspektasi frekuensi pada baris ke-i, kolom ke-j
A: Jenis Pekerjaan (i=1:Kantoran, i=2:Lapangan)
B: Tingkat Stres (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 :

\[ \log(\mu_{ij}) = \lambda \]

\[ + \lambda_1^A \, (\text{Kantoran}), \, \lambda_2^A \, (\text{Lapangan}) \]

\[ + \lambda_1^B \, (\text{Rendah}), \, \lambda_2^B \, (\text{Sedang}), \, \lambda_3^B \, (\text{Tinggi}) \]

\[ +\lambda_{ij}, (\text{interaksi jika ada}) \]

14.18 FITTING MODEL LOG-LINEAR DI R

# Membuat data frame dari tabel
tabel2x3 <- matrix(c(14, 18, 10, 12, 20, 18), nrow = 2, byrow = TRUE)
colnames(tabel2x3) <- c("Rendah", "Sedang", "Tinggi")
rownames(tabel2x3) <- c("kantoran", "Lapangan")
tabel2x3

##          Rendah Sedang Tinggi
## kantoran     14     18     10
## Lapangan     12     20     18

# Ubah menjadi data.frame untuk glm
data2x3 <- as.data.frame(as.table(tabel2x3))
colnames(data2x3) <- c("JenisPekerjaan", "TingkatStres", "Freq")
data2x3

##   JenisPekerjaan TingkatStres Freq
## 1       kantoran       Rendah   14
## 2       Lapangan       Rendah   12
## 3       kantoran       Sedang   18
## 4       Lapangan       Sedang   20
## 5       kantoran       Tinggi   10
## 6       Lapangan       Tinggi   18

# Model log-linear tanpa interaksi (asumsi independen)
fit_no_inter <- glm(Freq ~ JenisPekerjaan + TingkatStres, family = poisson, data = data2x3)
summary(fit_no_inter)

## 
## Call:
## glm(formula = Freq ~ JenisPekerjaan + TingkatStres, family = poisson, 
##     data = data2x3)
## 
## Coefficients:
##                        Estimate Std. Error z value Pr(>|z|)    
## (Intercept)             2.47398    0.22672  10.912   <2e-16 ***
## JenisPekerjaanLapangan  0.17435    0.20931   0.833    0.405    
## TingkatStresSedang      0.37949    0.25451   1.491    0.136    
## TingkatStresTinggi      0.07411    0.27235   0.272    0.786    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 5.1938  on 5  degrees of freedom
## Residual deviance: 1.8807  on 2  degrees of freedom
## AIC: 37.18
## 
## Number of Fisher Scoring iterations: 4

# Model log-linear dengan interaksi (untuk cek asosiasi)
fit_inter <- glm(Freq ~ JenisPekerjaan * TingkatStres, family = poisson, data = data2x3)
summary(fit_inter)

## 
## Call:
## glm(formula = Freq ~ JenisPekerjaan * TingkatStres, family = poisson, 
##     data = data2x3)
## 
## Coefficients:
##                                           Estimate Std. Error z value Pr(>|z|)
## (Intercept)                                 2.6391     0.2673   9.874   <2e-16
## JenisPekerjaanLapangan                     -0.1542     0.3934  -0.392    0.695
## TingkatStresSedang                          0.2513     0.3563   0.705    0.481
## TingkatStresTinggi                         -0.3365     0.4140  -0.813    0.416
## JenisPekerjaanLapangan:TingkatStresSedang   0.2595     0.5102   0.509    0.611
## JenisPekerjaanLapangan:TingkatStresTinggi   0.7419     0.5571   1.332    0.183
##                                              
## (Intercept)                               ***
## JenisPekerjaanLapangan                       
## TingkatStresSedang                           
## TingkatStresTinggi                           
## JenisPekerjaanLapangan:TingkatStresSedang    
## JenisPekerjaanLapangan:TingkatStresTinggi    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance:  5.1938e+00  on 5  degrees of freedom
## Residual deviance: -1.7764e-15  on 0  degrees of freedom
## AIC: 39.3
## 
## Number of Fisher Scoring iterations: 3

Interpretasi :

Untuk Model 1
- Tidak ada bukti yang cukup kuat untuk menyatakan ada hubungan antara Jenis Pekerjaan dan Tingkat Stres secara signifikan.
Untuk Model 2
- Interaksi untuk Lapangan × Tinggi sebesar 0.7419, menunjukkan bahwa pekerja lapangan dengan stres tinggi cenderung lebih banyak dibanding ekspektasi jika tidak ada interaksi. Namun tidak signifikan (p-value = 0.183).
Perbandingan Model
- Tidak terdapat bukti statistik yang kuat bahwa ada interaksi antara Jenis Pekerjaan dan Tingkat Stres. Model tanpa interaksi (asumsi independensi) cukup baik untuk menjelaskan hubungan antar variabel.

14.19 MODEL LOG LINEAR 3 ARAH

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

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

14.20 MODEL LOG-LINEAR UNTUK 3 ARAH

Model log-linear yang melibatkan tiga variabel kategorik (misal: X, Y, dan Z) dapat dibangun dalam berbagai bentuk model, tergantung pada tingkat interaksi yang ingin dimasukkan. Berikut adalah beberapa alternatif model log-linear yang umum digunakan:

1. Model Saturated \[ \log(\mu_{ijk}) = \lambda + \lambda^{X}_i + \lambda^{Y}_j + \lambda^{Z}_k + \lambda^{XY}_{ij} + \lambda^{XZ}_{ik} + \lambda^{YZ}_{jk} + \lambda^{XYZ}_{ijk} \] Model ini memuat semua kemungkinan interaksi, termasuk interaksi tiga arah (X, Y, dan Z).

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.

14.21 PENGUJIAN INTERAKSI DALAM MODEL LOG-LINEAR 3 ARAH

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

Pengujian Interaksi Tiga Arah (XYZ):
- Bandingkan model saturated dengan model homogenous.
Pengujian Interaksi Dua Arah (XY, XZ, YZ):
- Bandingkan model homogenous dengan model conditional.
- Bandingkan model conditional dengan model joint independence.
- Bandingkan model joint independence dengan model tanpa interaksi.

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

14.22 SOAL PRAKTIKUM

Tabel berikut menyajikan data hubungan antara Jenis Kelamin, Tingkat Pendidikan, dan Sikap terhadap Penggunaan Energi Nuklir sebagai sumber energi alternatif pada tahun 2023.

14.22.1 TABEL DATA SURVEI

Jenis Kelamin	Tingkat Pendidikan	Mendukung	Menolak	Total
Laki-laki	Rendah	40	20	60
Laki-laki	Menengah	55	25	80
Laki-laki	Tinggi	65	15	80
Perempuan	Rendah	30	25	55
Perempuan	Menengah	45	35	80
Perempuan	Tinggi	50	25	75

14.23 ANALISIS LOG-LINEAR UNTUK TABEL KONTINGENSI 3 ARAH

library("epitools")
library("DescTools")
library("lawstat")

## Warning: package 'lawstat' was built under R version 4.4.3

# Input data sesuai tabel praktikum
# Faktor-faktor
pendidikan <- factor(rep(c("Rendah", "Menengah", "Tinggi"), each = 4))
jenis_kelamin <- factor(rep(c("Laki-laki", "Perempuan"), each = 2, times = 3))
sikap <- factor(rep(c("Mendukung", "Menolak"), times = 6))

# Frekuensi sesuai tabel
frekuensi <- c(40, 20, 30, 25, 55, 25, 45, 35, 65, 15, 50, 25)

# Susun dalam data.frame
data_nuklir <- data.frame(
  Pendidikan = pendidikan,
  Jenis_Kelamin = jenis_kelamin,
  Sikap = sikap,
  Frekuensi = frekuensi
)

# Lihat datanya
data_nuklir

##    Pendidikan Jenis_Kelamin     Sikap Frekuensi
## 1      Rendah     Laki-laki Mendukung        40
## 2      Rendah     Laki-laki   Menolak        20
## 3      Rendah     Perempuan Mendukung        30
## 4      Rendah     Perempuan   Menolak        25
## 5    Menengah     Laki-laki Mendukung        55
## 6    Menengah     Laki-laki   Menolak        25
## 7    Menengah     Perempuan Mendukung        45
## 8    Menengah     Perempuan   Menolak        35
## 9      Tinggi     Laki-laki Mendukung        65
## 10     Tinggi     Laki-laki   Menolak        15
## 11     Tinggi     Perempuan Mendukung        50
## 12     Tinggi     Perempuan   Menolak        25

Membentuk Tabel Kontingensi 3 Arah

table3d <- xtabs(frekuensi ~ pendidikan + jenis_kelamin + sikap, data = data_nuklir)
ftable(table3d)

##                          sikap Mendukung Menolak
## pendidikan jenis_kelamin                        
## Menengah   Laki-laki                  55      25
##            Perempuan                  45      35
## Rendah     Laki-laki                  40      20
##            Perempuan                  30      25
## Tinggi     Laki-laki                  65      15
##            Perempuan                  50      25

Analisis Log-Linear: Tahap Pemodelan

Kita akan memodelkan tabel ini menggunakan beberapa model log-linear dan membandingkan kecocokan model (parsimonious model):

14.24 UJI MODEL INTERAKSI TIGA ARAH (SATURATED VS HOMOGENOUS)

14.24.1 PENENTUAN KATEGORI REFERENSI

##=============================##
# Penentuan kategori reference
##=============================##
jenis_kelamin <- relevel(jenis_kelamin, ref = "Perempuan")
sikap         <- relevel(sikap, ref = "Mendukung")
pendidikan    <- relevel(pendidikan, ref = "Tinggi")

Model Saturated Model log-linear saturated memasukkan semua interaksi hingga tiga arah :

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

# Model saturated
model_saturated <- glm(frekuensi ~ jenis_kelamin + sikap + pendidikan +
             jenis_kelamin*sikap + jenis_kelamin*pendidikan + sikap*pendidikan +
             jenis_kelamin*sikap*pendidikan,
             family = poisson(link = "log"))
summary(model_saturated)

## 
## Call:
## glm(formula = frekuensi ~ jenis_kelamin + sikap + pendidikan + 
##     jenis_kelamin * sikap + jenis_kelamin * pendidikan + sikap * 
##     pendidikan + jenis_kelamin * sikap * pendidikan, family = poisson(link = "log"))
## 
## Coefficients:
##                                                        Estimate Std. Error
## (Intercept)                                             3.91202    0.14142
## jenis_kelaminLaki-laki                                  0.26236    0.18811
## sikapMenolak                                           -0.69315    0.24495
## pendidikanMenengah                                     -0.10536    0.20548
## pendidikanRendah                                       -0.51083    0.23094
## jenis_kelaminLaki-laki:sikapMenolak                    -0.77319    0.37690
## jenis_kelaminLaki-laki:pendidikanMenengah              -0.06169    0.27530
## jenis_kelaminLaki-laki:pendidikanRendah                 0.02532    0.30613
## sikapMenolak:pendidikanMenengah                         0.44183    0.33286
## sikapMenolak:pendidikanRendah                           0.51083    0.36515
## jenis_kelaminLaki-laki:sikapMenolak:pendidikanMenengah  0.23605    0.50103
## jenis_kelaminLaki-laki:sikapMenolak:pendidikanRendah    0.26236    0.53887
##                                                        z value Pr(>|z|)    
## (Intercept)                                             27.662  < 2e-16 ***
## jenis_kelaminLaki-laki                                   1.395  0.16309    
## sikapMenolak                                            -2.830  0.00466 ** 
## pendidikanMenengah                                      -0.513  0.60812    
## pendidikanRendah                                        -2.212  0.02697 *  
## jenis_kelaminLaki-laki:sikapMenolak                     -2.051  0.04022 *  
## jenis_kelaminLaki-laki:pendidikanMenengah               -0.224  0.82268    
## jenis_kelaminLaki-laki:pendidikanRendah                  0.083  0.93409    
## sikapMenolak:pendidikanMenengah                          1.327  0.18438    
## sikapMenolak:pendidikanRendah                            1.399  0.16183    
## jenis_kelaminLaki-laki:sikapMenolak:pendidikanMenengah   0.471  0.63755    
## jenis_kelaminLaki-laki:sikapMenolak:pendidikanRendah     0.487  0.62635    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 7.1402e+01  on 11  degrees of freedom
## Residual deviance: 2.8422e-14  on  0  degrees of freedom
## AIC: 88.027
## 
## Number of Fisher Scoring iterations: 3

exp(model_saturated$coefficients)

##                                            (Intercept) 
##                                             50.0000000 
##                                 jenis_kelaminLaki-laki 
##                                              1.3000000 
##                                           sikapMenolak 
##                                              0.5000000 
##                                     pendidikanMenengah 
##                                              0.9000000 
##                                       pendidikanRendah 
##                                              0.6000000 
##                    jenis_kelaminLaki-laki:sikapMenolak 
##                                              0.4615385 
##              jenis_kelaminLaki-laki:pendidikanMenengah 
##                                              0.9401709 
##                jenis_kelaminLaki-laki:pendidikanRendah 
##                                              1.0256410 
##                        sikapMenolak:pendidikanMenengah 
##                                              1.5555556 
##                          sikapMenolak:pendidikanRendah 
##                                              1.6666667 
## jenis_kelaminLaki-laki:sikapMenolak:pendidikanMenengah 
##                                              1.2662338 
##   jenis_kelaminLaki-laki:sikapMenolak:pendidikanRendah 
##                                              1.3000000

Interpretasi :

Sikap terhadap energi nuklir tergantung pada jenis kelamin dan pendidikan, tapi pengaruh gabungan ketiganya (3-way interaction) tidak signifikan.
Ada indikasi laki-laki lebih kecil kemungkinan menolak energi nuklir dibanding perempuan, terutama pada kelompok tertentu. Pendidikan rendah juga dikaitkan dengan dukungan lebih tinggi terhadap energi nuklir.

14.24.2 RINGKASAN MODEL

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, sikap terhadap energi nuklir, dan tingkat pendidikan terhadap frekuensi responden.

14.24.3 HASIL ESTIMASI KOEFISIEN

Parameter	Estimate	Std. Error	z value	Pr(>z)
(Intercept)	3.92	0.15	27.67	<2e-16***
Laki-laki	0.27	0.19	1.4	0.164
Menolak	-0.7	0.25	-2.83	0.005**
Pendidikan Menengah	-0.11	0.21	-0.52	0.609
Pendidikan Rendah	-0.52	0.24	-2.22	0.027*
Laki-laki+Menolak	-0.78	0.38	-2.05	0.041*
Laki-laki+Menengah	-0.06	0.28	-0.23	0.823
Laki-laki+Rendah	0.03	0.31	0.09	0.935
Menolak+Menengah	0.45	0.34	1.33	0.185
Menolak+Rendah	0.52	0.37	1.4	0.162
Laki-laki+Menolak+Menengah	0.24	0.51	0.48	0.638
Laki-laki+Menolak+Rendah	0.27	0.54	0.49	0.627

14.24.4 INTERPRETASI KOEFISIEN

(Intercept): Rata-rata log jumlah kasus untuk kategori referensi (Perempuan, Menolak, Pendidikan Tinggi) adalah 3.92.
Menolak: Mereka yang menolak energi nuklir memiliki expected count sekitar 0.7 kali lipat lebih rendah dibanding yang mendukung (signifikan, p = 0.05).
Pendidikan Rendah: Pendidikan rendah dikaitkan dengan frekuensi yang lebih rendah.
Laki-laki & Menolak: Terdapat interaksi dua arah yang signifikan, dimana laki-laki cenderung lebih kecil menolak (signifikan, p = 0.041).
Interaksi dua & tiga arah: Interaksi tiga arah tidak signifikan sehingga kemungkinan model dapat disederhanakan.

14.24.5 GOODNESS-OF-FIT

Residual deviance ≈0 menandakan model saturated benar-benar fit terhadap data (seluruh variasi data dijelaskan oleh model).
AIC = 88.027 dapat digunakan untuk perbandingan dengan model yang lebih sederhana.

14.24.6 KESIMPULAN

Model saturated ini sangat fit dengan data, namun tidak semua parameter/interaksi signifikan.
Efek utama yang paling signifikan adalah:
- Menolak (expected count 0.7x lebih rendah dari yang mendukung)
- Pendidikan rendah (expected count 0.52x lebih rendah dari pendidikan tinggi)
- Tidak ditemukan bukti kuat interaksi tiga arah yang signifikan.
Tidak ditemukan bukti kuat interaksi dua atau tiga arah yang signifikan.
Model yang lebih sederhana (tanpa interaksi tiga arah) perlu dipertimbangkan untuk model final yang lebih parsimonious.

14.25 MODEL HOMOGENOUS

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

# Homogenous Model
model_homogenous <- glm(frekuensi ~ jenis_kelamin + sikap + pendidikan +
              jenis_kelamin*sikap + jenis_kelamin*pendidikan + sikap*pendidikan,
              family = poisson(link = "log"))
summary(model_homogenous)

## 
## Call:
## glm(formula = frekuensi ~ jenis_kelamin + sikap + pendidikan + 
##     jenis_kelamin * sikap + jenis_kelamin * pendidikan + sikap * 
##     pendidikan, family = poisson(link = "log"))
## 
## Coefficients:
##                                            Estimate Std. Error z value Pr(>|z|)
## (Intercept)                                3.936038   0.132812  29.636  < 2e-16
## jenis_kelaminLaki-laki                     0.219476   0.170504   1.287  0.19802
## sikapMenolak                              -0.766994   0.206914  -3.707  0.00021
## pendidikanMenengah                        -0.142769   0.186701  -0.765  0.44445
## pendidikanRendah                          -0.555607   0.208296  -2.667  0.00764
## jenis_kelaminLaki-laki:sikapMenolak       -0.602565   0.208410  -2.891  0.00384
## jenis_kelaminLaki-laki:pendidikanMenengah  0.005416   0.228917   0.024  0.98112
## jenis_kelaminLaki-laki:pendidikanRendah    0.104269   0.250471   0.416  0.67720
## sikapMenolak:pendidikanMenengah            0.546036   0.248057   2.201  0.02772
## sikapMenolak:pendidikanRendah              0.629803   0.267727   2.352  0.01865
##                                              
## (Intercept)                               ***
## jenis_kelaminLaki-laki                       
## sikapMenolak                              ***
## pendidikanMenengah                           
## pendidikanRendah                          ** 
## jenis_kelaminLaki-laki:sikapMenolak       ** 
## jenis_kelaminLaki-laki:pendidikanMenengah    
## jenis_kelaminLaki-laki:pendidikanRendah      
## sikapMenolak:pendidikanMenengah           *  
## sikapMenolak:pendidikanRendah             *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 71.40209  on 11  degrees of freedom
## Residual deviance:  0.30241  on  2  degrees of freedom
## AIC: 84.33
## 
## Number of Fisher Scoring iterations: 3

14.26 UJI HIPOTESIS : APAKAH ADA INTERAKSI 3 ARAH? (SATURATED VS HOMOGENOUS)

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

14.26.1 LANGKAH-LANGKAH PENGUJIAN

14.26.1.1 HIPOTESIS

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

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

14.26.1.2 HITUNG SELISIH DEVIANCE

# Deviance antar model
Deviance.model <- model_homogenous$deviance - model_saturated$deviance
Deviance.model

## [1] 0.3024145

14.26.1.3 HITUNG DERAJAT BEBAS

# Derajat bebas = db model homogenous - db model saturated
derajat.bebas <- (model_homogenous$df.residual - model_saturated$df.residual)
derajat.bebas

## [1] 2

14.26.1.4 CHI-SQUARE TABEL (α=0.05)

chi.tabel <- qchisq(1 - 0.05, df = derajat.bebas)
chi.tabel

## [1] 5.991465

14.26.1.5 KEPUTUSAN UJI

Keputusan <- ifelse(Deviance.model <= chi.tabel, "Terima H0", "Tolak H0")
Keputusan

## [1] "Terima H0"

Interpretasi : Pada taraf nyata 5%, belum cukup bukti untuk menolak H0 atau dapat dikatakan bahwa tidak ada interaksi tiga arah antara jenis kelamin, sikap terhadap energi nuklir, dan tingkat pendidikan.

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 = 0.302 - 0.00 = 0.302
- \(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 0.302 < 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, sikap terhadap energi nuklir, dan pendidikan.

14.27 UJI MODEL INTERAKSI DUA ARAH (HOMOGENOUS VS CONDITIONAL ON X)

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

# Conditional Association on X
model_conditional_X <- glm(frekuensi ~ jenis_kelamin + sikap + pendidikan +
              jenis_kelamin*sikap + jenis_kelamin*pendidikan,
              family = poisson(link = "log"))
summary(model_conditional_X)

## 
## Call:
## glm(formula = frekuensi ~ jenis_kelamin + sikap + pendidikan + 
##     jenis_kelamin * sikap + jenis_kelamin * pendidikan, family = poisson(link = "log"))
## 
## Coefficients:
##                                           Estimate Std. Error z value Pr(>|z|)
## (Intercept)                                3.79869    0.12873  29.509  < 2e-16
## jenis_kelaminLaki-laki                     0.26488    0.17543   1.510  0.13108
## sikapMenolak                              -0.38566    0.14059  -2.743  0.00608
## pendidikanMenengah                         0.06454    0.16073   0.402  0.68802
## pendidikanRendah                          -0.31015    0.17753  -1.747  0.08062
## jenis_kelaminLaki-laki:sikapMenolak       -0.59517    0.20659  -2.881  0.00397
## jenis_kelaminLaki-laki:pendidikanMenengah -0.06454    0.22546  -0.286  0.77469
## jenis_kelaminLaki-laki:pendidikanRendah    0.02247    0.24634   0.091  0.92731
##                                              
## (Intercept)                               ***
## jenis_kelaminLaki-laki                       
## sikapMenolak                              ** 
## pendidikanMenengah                           
## pendidikanRendah                          .  
## jenis_kelaminLaki-laki:sikapMenolak       ** 
## jenis_kelaminLaki-laki:pendidikanMenengah    
## jenis_kelaminLaki-laki:pendidikanRendah      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 71.4021  on 11  degrees of freedom
## Residual deviance:  7.3888  on  4  degrees of freedom
## AIC: 87.416
## 
## Number of Fisher Scoring iterations: 4

Pengujian Ada Tidaknya Interaksi Antara sikap dan pendidikan (Homogenous Model vs Conditional Association on X)

Hipotesis
- H0: \(\lambda_{jk}^{YZ}=0\) (Tidak ada interaksi antara sikap dan pendidikan)
- H0: \(\lambda_{jk}^{YZ}\neq0\) (Ada interaksi antara sikap dan pendidikan)
Tingkat Signifikansi

\(\alpha=5\%\)
Statistik Uji
- \(\Delta\)Deviance = Deviance model conditional on X - Deviance model homogenous = 7.086
- \(db\) = \(db\) model conditional on X - \(db\) model homogenous = 2
Daerah Penolakan

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

Karena 7.086 > 5.991, maka tolak H0

Interpretasi

Pada taraf nyata 5%, terdapat cukup bukti untuk menolak H0, atau dengan kata lain ada interaksi antara sikap dan pendidikan. Model terbaik yang terbentuk adalah model dengan menyertakan parameter interaksi (conditional terhadap jenis kelamin).

14.28 UJI MODEL INTERAKSI DUA ARAH (HOMOGENOUS VS CONDITIONAL ON Y)

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

# Conditional Association on Y
model_conditional_Y <- glm(frekuensi ~ jenis_kelamin + sikap + pendidikan +
              jenis_kelamin*sikap + sikap*pendidikan,
              family = poisson(link = "log"))
summary(model_conditional_Y)

## 
## Call:
## glm(formula = frekuensi ~ jenis_kelamin + sikap + pendidikan + 
##     jenis_kelamin * sikap + sikap * pendidikan, family = poisson(link = "log"))
## 
## Coefficients:
##                                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                           3.9208     0.1148  34.143  < 2e-16 ***
## jenis_kelaminLaki-laki                0.2469     0.1194   2.068 0.038643 *  
## sikapMenolak                         -0.7660     0.2075  -3.691 0.000223 ***
## pendidikanMenengah                   -0.1398     0.1367  -1.022 0.306705    
## pendidikanRendah                     -0.4964     0.1516  -3.275 0.001058 ** 
## jenis_kelaminLaki-laki:sikapMenolak  -0.5952     0.2066  -2.881 0.003966 ** 
## sikapMenolak:pendidikanMenengah       0.5452     0.2457   2.219 0.026474 *  
## sikapMenolak:pendidikanRendah         0.6142     0.2650   2.318 0.020440 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 71.40209  on 11  degrees of freedom
## Residual deviance:  0.51407  on  4  degrees of freedom
## AIC: 80.541
## 
## Number of Fisher Scoring iterations: 3

Pengujian Ada Tidaknya Interaksi Antara jenis kelamin dan pendidikan (Homogenous Model vs Conditional Association on Y)

Hipotesis
- H0: \(\lambda_{ik}^{XZ}=0\) (Tidak ada interaksi antara jenis kelamin dan pendidikan)
- H0: \(\lambda_{ik}^{XZ}\neq0\) (Ada interaksi antara jenis kelamin dan pendidikan)
Tingkat Signifikansi

\(\alpha=5\%\)
Statistik Uji
- \(\Delta\)Deviance = Deviance model conditional on Y - Deviance model homogenous = 0.212
- \(db\) = \(db\) model conditional on Y - \(db\) model homogenous = 2
Daerah Penolakan

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

Karena 0.212 < 5.991, maka terima H0
Interpretasi

Pada taraf nyata 5%, tidak terdapat cukup bukti untuk menolak H0, atau dengan kata lain tidak terdapat interaksi antara jenis kelamin dan pendidikan dan model homogenous lebih baik.

14.29 UJI MODEL INTERAKSI DUA ARAH (HOMOGENOUS VS CONDITIONAL ON Z)

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

# Conditional Association on Z
model_conditional_Z <- glm(frekuensi ~ jenis_kelamin + sikap + pendidikan +
              jenis_kelamin*pendidikan + sikap*pendidikan,
              family = poisson(link = "log"))
summary(model_conditional_Z)

## 
## Call:
## glm(formula = frekuensi ~ jenis_kelamin + sikap + pendidikan + 
##     jenis_kelamin * pendidikan + sikap * pendidikan, family = poisson(link = "log"))
## 
## Coefficients:
##                                           Estimate Std. Error z value Pr(>|z|)
## (Intercept)                                4.01900    0.12481  32.201  < 2e-16
## jenis_kelaminLaki-laki                     0.06454    0.16073   0.402   0.6880
## sikapMenolak                              -1.05605    0.18356  -5.753 8.76e-09
## pendidikanMenengah                        -0.10697    0.17840  -0.600   0.5488
## pendidikanRendah                          -0.50810    0.19837  -2.561   0.0104
## jenis_kelaminLaki-laki:pendidikanMenengah -0.06454    0.22546  -0.286   0.7747
## jenis_kelaminLaki-laki:pendidikanRendah    0.02247    0.24634   0.091   0.9273
## sikapMenolak:pendidikanMenengah            0.54523    0.24569   2.219   0.0265
## sikapMenolak:pendidikanRendah              0.61422    0.26496   2.318   0.0204
##                                              
## (Intercept)                               ***
## jenis_kelaminLaki-laki                       
## sikapMenolak                              ***
## pendidikanMenengah                           
## pendidikanRendah                          *  
## jenis_kelaminLaki-laki:pendidikanMenengah    
## jenis_kelaminLaki-laki:pendidikanRendah      
## sikapMenolak:pendidikanMenengah           *  
## sikapMenolak:pendidikanRendah             *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 71.4021  on 11  degrees of freedom
## Residual deviance:  8.7764  on  3  degrees of freedom
## AIC: 90.804
## 
## Number of Fisher Scoring iterations: 4

Pengujian Ada Tidaknya Interaksi Antara jenis kelamin dan sikap (Homogenous Model vs Conditional Association on Z)

Hipotesis
- H0: \(\lambda_{ij}^{XY}=0\) (Tidak ada interaksi antara jenis kelamin dan sikap)
- H0: \(\lambda_{ij}^{XY}\neq0\) (Ada interaksi antara jenis kelamin dan sikap)
Tingkat Signifikansi

\(\alpha=5\%\)
Statistik Uji
- \(\Delta\)Deviance = Deviance model conditional on Y - Deviance model homogenous = 8.776
- \(db\) = \(db\) model homogenous - \(db\) model saturated = 1
Daerah Penolakan

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

Karena 8.776 > 3.841, maka tolak H0
Interpretasi

Dengan taraf nyata 5%, menolak H0 atau dapat dikatakan bahwa ada interaksi antara jenis kelamin dan pendapat tentang sikap. Dengan kata lain, model yang terbentuk adalah model parameter interaksi (conditional terhadap pendidikan).

14.30 PEMILIHAN MODEL TERBAIK

14.31 RINGKASAN MODEL LOG-LINEAR

Model	Deviance	df	AIC
Saturated	0.000	0	88.027
Homogenous	0.302	2	84.33
Conditional on X	7.388	4	87.416
Conditional on Y	0.514	4	80.541
Conditional on Z	8.776	3	90.804

14.32 RINGKASAN PENGUJIAN 2 ARAH DAN 3 ARAH

Interaksi	Pengujian	\(\Delta\)Deviance	\(\Delta\)df	Chi-Square Tabel	Keputusan	Keterangan
XYZ	Saturated vs Homogenous	0.302	2	5.991	Terima H0	tidak ada interaksi
YZ	Conditional on X vs Homogenous	7.086	2	5.991	Tolak H0	ada interaksi
XZ	Conditional on Y vs Homogenous	0.212	2	5.991	Terima H0	tidak ada interaksi
XY	Conditional on Z vs Homogenous	8,776	1	3.841	Tolak H0	ada interaksi

14.33 KESIMPULAN PEMILIHAN MODEL TERBAIK

Dari hasil di atas diketahui bahwa model terbaik yang diperoleh adalah model homogenous adalah:

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

Model terbaik adalah model log-linear tanpa interaksi tiga arah dan hanya semua interaksi dua arah.

14.34 MODEL TERBAIK

Model terbaik dipilih berdasarkan pengujian interaksi yang signifikan, yaitu model homogenous:

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

# Model Terbaik
bestmodel <- glm(frekuensi ~ jenis_kelamin + sikap + pendidikan +
              jenis_kelamin*sikap + jenis_kelamin*pendidikan + sikap*pendidikan,
              family = poisson(link = "log"))
summary(bestmodel)

## 
## Call:
## glm(formula = frekuensi ~ jenis_kelamin + sikap + pendidikan + 
##     jenis_kelamin * sikap + jenis_kelamin * pendidikan + sikap * 
##     pendidikan, family = poisson(link = "log"))
## 
## Coefficients:
##                                            Estimate Std. Error z value Pr(>|z|)
## (Intercept)                                3.936038   0.132812  29.636  < 2e-16
## jenis_kelaminLaki-laki                     0.219476   0.170504   1.287  0.19802
## sikapMenolak                              -0.766994   0.206914  -3.707  0.00021
## pendidikanMenengah                        -0.142769   0.186701  -0.765  0.44445
## pendidikanRendah                          -0.555607   0.208296  -2.667  0.00764
## jenis_kelaminLaki-laki:sikapMenolak       -0.602565   0.208410  -2.891  0.00384
## jenis_kelaminLaki-laki:pendidikanMenengah  0.005416   0.228917   0.024  0.98112
## jenis_kelaminLaki-laki:pendidikanRendah    0.104269   0.250471   0.416  0.67720
## sikapMenolak:pendidikanMenengah            0.546036   0.248057   2.201  0.02772
## sikapMenolak:pendidikanRendah              0.629803   0.267727   2.352  0.01865
##                                              
## (Intercept)                               ***
## jenis_kelaminLaki-laki                       
## sikapMenolak                              ***
## pendidikanMenengah                           
## pendidikanRendah                          ** 
## jenis_kelaminLaki-laki:sikapMenolak       ** 
## jenis_kelaminLaki-laki:pendidikanMenengah    
## jenis_kelaminLaki-laki:pendidikanRendah      
## sikapMenolak:pendidikanMenengah           *  
## sikapMenolak:pendidikanRendah             *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 71.40209  on 11  degrees of freedom
## Residual deviance:  0.30241  on  2  degrees of freedom
## AIC: 84.33
## 
## Number of Fisher Scoring iterations: 3

Model	AIC
Saturated	88.027
Homogenous	84.33
Conditional on X	87.416
Conditional on Y	80.541
Conditional on Z	90.804
Model Terbaik	84.33

Dari hasil di atas, terlihat bahwa model terbaik memiliki AIC yang rendah dibandingkan dengan model lainnya seperti saturated dan conditional model. Meskipun conditional on Y memiliki AIC yang lebih kecil, tetapi hasil uji menunjukkan bahwa model homogenous lebih baik dibandingkan ketika interaksi dua arah ditambahkan.

14.35 INTERPRETASI KOEFISIEN MODEL TERBAIK

# Interpretasi koefisien model terbaik
data.frame(
  koef = bestmodel$coefficients,
  exp_koef = exp(bestmodel$coefficients)
)

##                                                   koef   exp_koef
## (Intercept)                                3.936037898 51.2152786
## jenis_kelaminLaki-laki                     0.219475787  1.2454237
## sikapMenolak                              -0.766994483  0.4644068
## pendidikanMenengah                        -0.142769458  0.8669539
## pendidikanRendah                          -0.555606817  0.5737240
## jenis_kelaminLaki-laki:sikapMenolak       -0.602565279  0.5474056
## jenis_kelaminLaki-laki:pendidikanMenengah  0.005415838  1.0054305
## jenis_kelaminLaki-laki:pendidikanRendah    0.104269119  1.1098991
## sikapMenolak:pendidikanMenengah            0.546035546  1.7263952
## sikapMenolak:pendidikanRendah              0.629802612  1.8772400

14.36 INTERPRETASI KOEFISIEN MODEL TERBAIK

\(\exp(0.219) = 1.245 \rightarrow nilai odds\)

Tanpa memperhatikan faktor lain, peluang seseorang berjenis kelamin laki-laki adalah 1.245 kali dibandingkan perempuan.
\(\exp(-0.767) = 0.464 \rightarrow nilai odds\)

Tanpa memperhatikan faktor lain, peluang seseorang menolak energi nuklir adalah 0.464 kali dibandingkan yang menerima.
\(\exp(-0.143) = 0.867 \rightarrow nilai odds\)

Tanpa memperhatikan faktor lain, peluang seseorang berpendidikan menengah adalah 0.867 kali dibandingkan yang tinggi.
\(\exp(-0.556) = 0.574 \rightarrow nilai odds\)

Tanpa memperhatikan faktor lain, peluang seseorang berpendidikan rendah adalah 0.574 kali dibandingkan yang tinggi.

14.37 NILAI DUGAAN MODEL TERBAIK

# Fitted values dari model terbaik
data.frame(
  Pendidikan = pendidikan,
  Jenis_Kelamin = jenis_kelamin,
  Sikap = sikap,
  fitted = bestmodel$fitted.values
)

##    Pendidikan Jenis_Kelamin     Sikap   fitted
## 1      Rendah     Laki-laki Mendukung 40.61656
## 2      Rendah     Laki-laki   Menolak 19.38344
## 3      Rendah     Perempuan Mendukung 29.38344
## 4      Rendah     Perempuan   Menolak 25.61656
## 5    Menengah     Laki-laki Mendukung 55.59871
## 6    Menengah     Laki-laki   Menolak 24.40129
## 7    Menengah     Perempuan Mendukung 44.40129
## 8    Menengah     Perempuan   Menolak 35.59871
## 9      Tinggi     Laki-laki Mendukung 63.78472
## 10     Tinggi     Laki-laki   Menolak 16.21528
## 11     Tinggi     Perempuan Mendukung 51.21528
## 12     Tinggi     Perempuan   Menolak 23.78472

15 REFERENSI

Agresti, A. (2018). An Introduction to Categorical Data Analysis. John Wiley & Sons.
Hosmer, D. W., Lemeshow, S., & Sturdivant, R. X. (2013). Applied Logistic Regression. John Wiley & Sons.
Wibisono, E. (2019). Statistik Terapan untuk Ilmu Sosial. Penerbit Andi.
Field, A. (2017). Discovering Statistics Using IBM SPSS Statistics. SAGE Publications.
Sudjana, T. (2005). Metode Statistika. Tarsito.
James, G., Witten, D., Hastie, T., & Tibshirani, R. (2013). An Introduction to Statistical Learning. Springer.
Han, J., Kamber, M., & Pei, J. (2011). Data Mining: Concepts and Techniques. Morgan Kaufmann.
Breiman, L. (2001). Random Forests. Machine Learning, 45(1), 5-32.
Walpole, R. E., Myers, R. H., Myers, S. L., & Ye, K. (2012). Probability and Statistics for Engineers and Scientists (9th ed.). Pearson.
Ross, S. M. (2014). Introduction to Probability and Statistics for Engineers and Scientists (5th ed.). Academic Press.
Verzani, J. (2005). Using R for Introductory Statistics. CRC Press.
R Documentation: https://stat.ethz.ch/R-manual/
Song, J. W., & Chung, K. C. (2010). “Observational studies: Cohort and case-control studies.” Plastic and reconstructive surgery, 126(6), 2234–2242. DOI: 10.1097/PRS.0b013e3181f44abc
Mann, C. J. (2003). “Observational research methods. Research design II: cohort, cross sectional, and case-control studies.” Emergency medicine journal, 20(1), 54–60. DOI: 10.1136/emj.20.1.54
Agresti, A. (2018). Statistical Methods for the Social Sciences. Pearson.
McHugh, M. L. (2013). The chi-square test of independence. Biochemia Medica, 23(2), 143–149.
Altman, D. G. (1991). Practical Statistics for Medical Research. Chapman & Hall.
Kirkwood, B. R., & Sterne, J. A. C. (2003). Essential Medical Statistics. Blackwell Science.
Zhang, J., & Yu, K. F. (1998). What’s the relative risk? JAMA, 280(19), 1690–1691. https://doi.org/10.1001/jama.280.19.1690
Fleiss, J. L., Levin, B., & Paik, M. C. (2003). Statistical Methods for Rates and Proportions. Wiley-Interscience.
Agresti, A. (2002). Categorical Data Analysis (2nd ed.). Wiley.
Walpole, R. E., Myers, R. H., Myers, S. L., & Ye, K. (2012). Probability and Statistics for Engineers and Scientists. Pearson.
Agresti, A. (2018). Statistical Methods for the Social Sciences. Pearson.
Newbold, P., Carlson, W. L., & Thorne, B. (2012). Statistics for Business and Economics.
Agresti, A. (2013). Categorical Data Analysis. Wiley.
Mehta, C. R., & Patel, N. R. (1983). A network algorithm for performing Fisher’s exact test in r×c contingency tables. Journal of the American Statistical Association.
Friendly, M. (2000). Visualizing Categorical Data. SAS Institute Inc.
Sharpe, D. (2015). Chi-Square Test is Statistically Significant: Now What?. Practical Assessment, Research, and Evaluation.
Kim, H. Y. (2017). Statistical notes for clinical researchers: Chi-squared test and Fisher’s exact test. Restorative Dentistry & Endodontics, 42(2), 152–155. https://doi.org/10.5395/rde.2017.42.2.152
Freeman, J. V., et al. (2012). Statistical methods for categorical data analysis. CRC Press.
Fisher, R. A. (1922). On the Interpretation of χ² from Contingency Tables, and the Calculation of P. Journal of the Royal Statistical Society, 85(1), 87-94.
Haberman, S. J. (1973). “The analysis of residuals in cross-classified tables.” Biometrics, 29(1), 205–220.
Friendly, M. (1994). “Mosaic displays for multi-way contingency tables.” Journal of the American Statistical Association, 89(425), 190–200.
Agresti, A. (2007). An Introduction to Categorical Data Analysis. John Wiley & Sons.
Wickens, T. D. (1989). Multiway Contingency Table Analysis for the Social Sciences. Erlbaum.
Breslow, N. E., & Day, N. E. (1980). Statistical Methods in Cancer Research. WHO/IARC.
Wagner, C. H. (1982). Simpson’s paradox in real life. The American Statistician.
Pearl, J. (2014). Understanding Simpson’s Paradox. The American Statistician.
Goodman, L. A. (1970). The multivariate analysis of three-way contingency tables. Journal of the Royal Statistical Society, Series B.
Lauritzen, S. L. (1996). Graphical Models. Oxford University Press.
Pearl, J. (2009). Causality: Models, Reasoning and Inference. Cambridge University Press.
Dawid, A. P. (1979). “Conditional Independence in Statistical Theory”. Journal of the Royal Statistical Society, Series B.
Pearl, J. (1988). Probabilistic Reasoning in Intelligent Systems: Networks of Plausible Inference. San Francisco: Morgan Kaufmann Publishers.
Spirtes, P., Glymour, C., & Scheines, R. (2000). Causation, Prediction, and Search. Cambridge, MA: MIT Press.
Rayner, J. C. W., & Best, D. J. (2018). Extensions to the Cochran–Mantel–Haenszel mean scores and correlation tests. Computational Statistics & Data Analysis, 120, 211-222. https://doi.org/10.1016/j.csda.2018.03.015
NCSS. (n.d.). Tests for two proportions in a stratified design (Cochran-Mantel-Haenszel test). Retrieved from https://www.ncss.com/wp-content/themes/ncss/pdf/Procedures/PASS/Tests_for_Two_Proportions_in_a_Stratified_Design-Cochran-Mantel-Haenszel_Test.pdf
J. T., & Wallenstein, S. (1988). The power of the Mantel-Haenszel test. Statistics in Medicine, 7(6), 757-764. https://doi.org/10.1002/sim.4780070610
Breslow, N. E., & Day, N. E. (1980). Statistical methods in cancer research: Volume I - The analysis of case-control studies. International Agency for Research on Cancer.
Fisher, R. A. (1935). The Design of Experiments. Oliver and Boyd.
Mehta, C. R., & Patel, N. R. (1983). A Network Algorithm for Performing Fisher’s Exact Test in r × c Contingency Tables. Journal of the American Statistical Association, 78(382), 427–434.
Howell, D. C. (2012). Statistical Methods for Psychology (8th ed.). Cengage Learning.

Analisis Data Kategori

Fazila Azra Anggina

1 PENDAHULUAN

1.1 TUJUAN

1.2 JENIS

1.3 PERBEDAAN DENGAN DATA KUANTITATIF

1.4 IMPLEMENTASI

2 METODE

2.1 REGRESI LOGISTIK

2.2 ANALISIS KORESPONDENSI

2.3 RANDOM FOREST

3 DISTRIBUSI PROBABILITAS

3.1 DISTRIBUSI BERNOULLI

3.2 DISTRIBUSI BINOMIAL

3.3 DISTRIBUSI MULTINOMIAL

3.4 DISTRIBUSI POISSON

4 DESAIN SAMPLING

4.1 PROSPEKTIF SAMPLING

4.1.1 EKSPERIMEN

4.1.2 STUDI KOHORT PROSPEKTIF

4.2 RETROSPEKTIF SAMPLING

4.2.1 STUDI KASUS-KONTROL

4.2.2 STUDI KOHORT RETROSPEKTIF

5 TABEL KONTINGENSI 2X2

5.1 DISTRIBUSI FREKUENSI DAN PROPORSI

5.1.1 FREKUENSI ABSOLUT

5.1.2 FREKUENSI RELATIF

5.1.3 DISTRIBUSI MARGINAL

5.1.4 DISTRIBUSI KONDISIONAL

5.2 DISTRIBUSI PELUANG DALAM TABEL KONTINGENSI 2X2

5.2.1 PELUANG BERSAMA (JOINT PROBABILITY)

5.2.2 PELUANG MARGINAL (MARGINAL PROBABILITY)

5.2.3 PELUANG BERSYARAT (CONDITIONAL PROBABILITY)

5.3 UKURAN ASOSIASI DALAM TABEL KONTINGENSI 2X2

5.3.1 RISK DIFFERENCE (PERBEDAAN RISIKO)

5.3.2 RELATIVE RISK (RISIKO RELATIF)

5.3.3 ODDS RATIO (RASIO ODDS)

6 INFERENSI TABEL KONTINGENSI 2X2

6.1 ESTIMASI

6.1.1 ESTIMASI TITIK

6.1.2 ESTIMASI INTERVAL

6.2 UJI HIPOTESIS

6.2.1 UJI PROPORSI

6.2.2 UJI ASOSIASI

6.2.3 UJI INDEPENDENSI

6.2.3.1 UJI CHI-SQUARE

6.2.3.2 UJI FISHER’S EXACT

6.3 ANALISIS RESIDU

6.3.1 JENIS RESIDUAL DAN DETEKSI OUTLIER

6.3.1.1 PEARSON RESIDUAL (RESIDU STANDAR)

6.3.1.2 STANDARDIZED RESIDUAL (Z-SCORE RESIDUAL)

6.3.1.3 DETEKSI OUTLIER

7 TABEL KONTINGENSI TIGA ARAH

7.1 TABEL PARSIAL DAN MARGINAL

7.2 DISTRIBUSI PELUANG

7.2.1 PELUANG BERSAMA (JOINT PROBABILITY)

7.2.2 PELUANG MARGINAL (MARGINAL PROBABILITY)

7.2.3 PELUANG BERSYARAT (CONDITIONAL PROBABILITY)

7.3 UKURAN ASOSIASI DALAM TABEL KONTINGENSI 3X3

7.3.1 RISK DIFFERENCE (PERBEDAAN RISIKO)

7.3.2 RELATIVE RISK (RISIKO RELATIF)

7.3.3 ODDS RATIO (RASIO ODDS)

7.4 INDEPENDENSI BERSYARAT

7.5 INFERENSI TABEL KONTINGENSI TIGA ARAH

7.5.1 COCHRAN-MANTEL-HAENSZEL

7.5.2 ODDS RATIO BERSAMA

7.5.3 UJI HOMOGENITAS ODDS RATIO DENGAN BRESLOW-DAY

8 GENERALIZED LINEAR MODEL (GLM)

8.1 MODEL REGRESI LOGISTIK

8.1.1 ESTIMASI PARAMETER

8.2 MODEL REGRESI POISSON

8.2.1 ESTIMASI PARAMETER

9 INFERENSI GENERALIZED LINEAR MODEL (GLM)

9.1 EKSPEKTASI DAN VARIANS

9.1.1 ESTIMASI PARAMETER

9.2 DIAGNOSTIK MODEL GLM

9.3 DETAIL METODE ESTIMASI DAN INFERENSI REGRESI LOGISTIK

9.3.1 ESTIMASI DENGAN NEWTON-RAPHSON

9.4 UJI WALD

9.5 UJI LIKELIHOOD RATIO