Ebook Analisi Data Kategori

1. Pendahuluan - Analisis Data Kategorik pada Kasus Insomnia

Analisis data kategorik merupakan salah satu pendekatan penting dalam ilmu statistik yang digunakan untuk memahami hubungan antar variabel yang bersifat kategorik. Dalam konteks penelitian sosial maupun medis, sering kali data yang dikumpulkan bukan dalam bentuk angka kontinu, melainkan berupa kategori. Contohnya adalah variabel seperti jenis kelamin, status merokok, kebiasaan mengonsumsi kopi, atau kondisi kesehatan seperti insomnia.

Pada bab ini, akan diperkenalkan latar belakang analisis data kategorik melalui studi kasus yang relevan, yaitu hubungan antara konsumsi kopi dan kejadian insomnia pada individu dengan mempertimbangkan faktor jenis kelamin sebagai variabel pengganggu potensial.

1.1 Latar Belakang Masalah

Insomnia merupakan salah satu gangguan tidur yang banyak dikeluhkan oleh masyarakat modern. Salah satu faktor yang sering diasosiasikan dengan insomnia adalah konsumsi kafein, khususnya dalam bentuk kopi. Banyak orang percaya bahwa minum kopi dapat menyebabkan kesulitan tidur. Namun, apakah benar bahwa semua orang yang minum kopi akan mengalami insomnia? Apakah pengaruh ini sama untuk laki-laki dan perempuan?

Analisis data kategorik memberikan kita alat untuk menjawab pertanyaan-pertanyaan ini secara sistematis dan berbasis data.

1.2 Tujuan Penelitian

Tujuan utama dari analisis ini adalah untuk:

Tujuan utama dari analisis ini adalah untuk mengetahui apakah terdapat hubungan antara konsumsi kopi dan kejadian insomnia. Hal ini penting untuk menguji anggapan umum bahwa kopi dapat memicu gangguan tidur.

-Selanjutnya, analisis ini ingin melihat apakah pengaruh tersebut konsisten pada laki-laki dan perempuan. Ini membantu mengetahui apakah jenis kelamin memengaruhi pola hubungan antara kopi dan insomnia.

Terakhir, tujuan lainnya adalah memahami apakah hubungan antara kopi dan insomnia berubah ketika mempertimbangkan variabel ketiga, yaitu gender. Ini penting untuk mengidentifikasi adanya efek perancu atau interaksi antar variabel.

1.3 Pentingnya Analisis Kategorik

Analisis data kategorik penting karena banyak variabel dalam bidang sosial dan kesehatan bersifat nominal atau ordinal, bukan numerik. Metode seperti regresi linear tidak sesuai karena tidak memenuhi asumsi dasar untuk data kategorik. Oleh karena itu, diperlukan pendekatan khusus seperti tabel kontingensi, regresi logistik, dan model log-linear.

Teknik-teknik ini memungkinkan kita mengevaluasi hubungan antar variabel kategorik, mengontrol efek dari variabel lain, serta mengidentifikasi pola interaksi. Dengan demikian, analisis kategorik memberikan landasan statistik yang kuat dalam menjawab pertanyaan-pertanyaan penelitian yang melibatkan data diskrit.

Melalui pendekatan ini, kita tidak hanya dapat melihat hubungan antara dua variabel, tetapi juga mengontrol efek variabel lain dan menguji interaksi antar variabel secara menyeluruh.

1.4 Struktur Buku

Struktur buku ini disusun secara bertahap, dimulai dari analisis dua arah menggunakan tabel kontingensi sederhana. Selanjutnya, pembahasan akan berkembang ke analisis tiga arah, di mana variabel ketiga turut diperhitungkan dalam hubungan antar dua variabel utama.

Kemudian, pembaca akan diperkenalkan pada pendekatan model yang lebih kompleks, seperti Generalized Linear Model (GLM), untuk mengakomodasi kompleksitas hubungan antar variabel. Setiap bagian didukung dengan simulasi data, analisis statistik, dan interpretasi hasil.

1.5 Contoh Simulasi Data Awal

Sebagai pengantar, data simulasi disiapkan untuk merepresentasikan tiga variabel utama: konsumsi kopi, kejadian insomnia, dan jenis kelamin. Data ini dibuat secara acak untuk menggambarkan kemungkinan distribusi pada populasi nyata.

set.seed(1)
df_intro <- data.frame(
  Kopi = sample(c("Ya", "Tidak"), 100, replace = TRUE),
  Insomnia = sample(c("Ya", "Tidak"), 100, replace = TRUE),
  Gender = sample(c("Laki", "Perempuan"), 100, replace = TRUE)
)
head(df_intro)

##    Kopi Insomnia    Gender
## 1    Ya    Tidak Perempuan
## 2 Tidak    Tidak Perempuan
## 3    Ya    Tidak      Laki
## 4    Ya       Ya Perempuan
## 5 Tidak       Ya Perempuan
## 6    Ya       Ya Perempuan

1.6 Visualisasi Awal Data

Visualisasi awal digunakan untuk melihat gambaran umum hubungan antara konsumsi kopi, insomnia, dan gender. Grafik batang proporsi memberikan ilustrasi perbandingan kejadian insomnia berdasarkan kategori kopi dan dibedakan menurut jenis kelamin.

Dari visualisasi ini, kita bisa mengamati pola awal yang akan menjadi dasar eksplorasi lebih lanjut dalam analisis statistik. Ini membantu memunculkan dugaan awal atau hipotesis mengenai hubungan antar variabel.

library(ggplot2)
ggplot(df_intro, aes(x = Kopi, fill = Insomnia)) +
  geom_bar(position = "fill") +
  facet_wrap(~Gender) +
  labs(y = "Proporsi", title = "Proporsi Insomnia Berdasarkan Konsumsi Kopi dan Gender")

Bab ini memberikan fondasi yang kuat untuk memahami konteks dan motivasi analisis. Dengan memahami struktur data dan tujuan dari analisis kategorik, kita dapat lebih siap memasuki bab-bab selanjutnya yang akan membahas metode dan aplikasi secara lebih teknis dan mendalam.

2. Metode dalam Analisis Data Kategori

Bab ini membahas empat metode utama dalam analisis data kategori: tabel kontingensi, regresi logistik, correspondence analysis (CA), dan decision tree. Setiap metode akan dilengkapi dengan rumus, simulasi, dan interpretasi. Kita gunakan contoh kasus: apakah konsumsi kopi berpengaruh terhadap insomnia.

2.1 Tabel Kontingensi dan Uji Chi-Square

Tabel kontingensi digunakan untuk menyajikan data kategorik dalam bentuk matriks frekuensi berdasarkan dua variabel. Dalam kasus ini, tabel menunjukkan distribusi jumlah orang yang mengalami atau tidak mengalami insomnia, berdasarkan status konsumsi kopi.

Setelah tabel dibentuk, uji Chi-Square digunakan untuk menguji apakah ada hubungan yang signifikan antara kedua variabel. Uji ini membandingkan nilai frekuensi yang diamati dengan nilai yang diharapkan jika tidak ada hubungan (independen).

Contoh Data (Hipotetik)

Konsumsi Kopi	Insomnia	Tidak Insomnia
Ya	45	15
Tidak	30	40

kopi_insomnia <- matrix(c(45, 15, 30, 40), nrow = 2, byrow = TRUE)
dimnames(kopi_insomnia) <- list(
  "Konsumsi Kopi" = c("Ya", "Tidak"),
  "Insomnia" = c("Ya", "Tidak")
)
addmargins(kopi_insomnia)

##              Insomnia
## Konsumsi Kopi Ya Tidak Sum
##         Ya    45    15  60
##         Tidak 30    40  70
##         Sum   75    55 130

Rumus Chi-Square

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

Dengan: - $O_{ij}$: nilai observasi - $E_{ij} = \frac{(\text{jumlah baris}) \times (\text{jumlah kolom})}{\text{total}}$: nilai ekspektasi

observed <- kopi_insomnia
expected <- outer(rowSums(observed), colSums(observed)) / sum(observed)
chisq_manual <- sum((observed - expected)^2 / expected)
chisq_manual

## [1] 13.67532

Uji Statistik di R

chisq.test(kopi_insomnia, correct = FALSE)

## 
##  Pearson's Chi-squared test
## 
## data:  kopi_insomnia
## X-squared = 13.675, df = 1, p-value = 0.0002173

Interpretasi

Jika nilai $p < 0.05$, maka kita tolak $H_0$ dan menyimpulkan bahwa ada hubungan antara konsumsi kopi dan insomnia.

2.2 Regresi Logistik

Regresi logistik adalah metode yang digunakan untuk memodelkan hubungan antara variabel prediktor (misalnya konsumsi kopi) dan variabel respons biner (misalnya status insomnia: ya/tidak). Model ini menghitung peluang seseorang mengalami insomnia berdasarkan status konsumsi kopi.

Rumus Model

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

set.seed(123)
konsumsi <- c(rep(1, 60), rep(0, 70)) # 60 peminum kopi, 70 tidak
insomnia <- c(rbinom(60, 1, 0.75), rbinom(70, 1, 0.43))
data <- data.frame(insomnia, konsumsi)
model <- glm(insomnia ~ konsumsi, data = data, family = binomial)
summary(model)

## 
## Call:
## glm(formula = insomnia ~ konsumsi, family = binomial, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)   
## (Intercept)  -0.2877     0.2415  -1.191  0.23361   
## konsumsi      1.2157     0.3747   3.244  0.00118 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 178.24  on 129  degrees of freedom
## Residual deviance: 167.14  on 128  degrees of freedom
## AIC: 171.14
## 
## Number of Fisher Scoring iterations: 4

Odds Ratio

exp(coef(model))

## (Intercept)    konsumsi 
##    0.750000    3.372549

Interpretasi

Koefisien positif menunjukkan bahwa konsumsi kopi meningkatkan peluang mengalami insomnia.
Odds ratio > 1 berarti kemungkinan insomnia lebih tinggi pada peminum kopi.

2.3 Correspondence Analysis (CA)

Correspondence Analysis (CA) adalah metode eksploratif yang digunakan untuk menganalisis hubungan antara dua variabel kategorik dalam bentuk visual. Teknik ini memungkinkan kita untuk menyajikan baris dan kolom dari tabel kontingensi dalam ruang berdimensi rendah, biasanya dua dimensi.

Dalam konteks kasus ini, CA membantu memetakan hubungan antara konsumsi kopi dan insomnia dengan melihat pola kedekatan antar kategori. Semakin dekat dua kategori pada grafik hasil CA, semakin besar asosiasi di antara keduanya.

2.4 Pohon Keputusan (Decision Tree)

Pohon keputusan adalah metode sederhana untuk membuat prediksi berdasarkan data. Dalam contoh ini, kita menggunakannya untuk memprediksi apakah seseorang mengalami insomnia atau tidak, dengan melihat apakah orang tersebut minum kopi.

Prosesnya bekerja seperti menjawab pertanyaan “ya” atau “tidak” secara bertahap. Setiap percabangan menunjukkan kondisi tertentu, misalnya “Apakah minum kopi?”, dan setiap ujung cabang (daun) menunjukkan hasil akhir, seperti “kemungkinan insomnia tinggi” atau “tidak insomnia”

library(rpart)
library(rpart.plot)

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

tree_model <- rpart(as.factor(insomnia) ~ konsumsi, data = data, method = "class")
rpart.plot(tree_model)

Secara keseluruhan, keempat metode dalam bab ini saling melengkapi dalam menganalisis data kategorik dan memberikan hasil yang konsisten: ada indikasi bahwa konsumsi kopi memang berkaitan dengan gangguan tidur.

3. Distribusi Probabilitas dalam Data Kategori

Dalam konteks konsumsi kopi dan insomnia, kita dapat memahami fenomena tersebut melalui beberapa distribusi probabilitas diskret, yaitu: Bernoulli, Binomial, Multinomial, dan Poisson. Masing-masing akan digunakan untuk menjelaskan tipe data kategorik dan hubungannya dengan kejadian insomnia.

3.1 Distribusi Bernoulli

Distribusi Bernoulli digunakan ketika hasil yang diamati hanya dua kemungkinan, misalnya: seseorang mengalami insomnia (1) atau tidak (0). Ini cocok untuk menggambarkan kejadian yang bersifat ya atau tidak.

Dalam kasus ini, kita bisa mengasumsikan bahwa peluang seseorang mengalami insomnia adalah 70%. Dengan menggunakan distribusi Bernoulli, kita bisa mensimulasikan data untuk melihat bagaimana hasil insomnia akan muncul jika peluangnya seperti itu.

Rumus

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

set.seed(1)
bernoulli_samples <- rbinom(20, size = 1, prob = 0.7)
mean(bernoulli_samples)  # Estimasi peluang insomnia

## [1] 0.65

Interpretasi

Jika peluang insomnia pada peminum kopi diasumsikan 70%, maka simulasi Bernoulli ini meniru keputusan insomnia pada 20 peminum secara acak.

3.2 Distribusi Binomial

Distribusi Binomial digunakan untuk menghitung jumlah kejadian tertentu dalam sejumlah percobaan tetap. Misalnya, dari 60 orang peminum kopi, berapa banyak yang mungkin mengalami insomnia jika peluangnya 70% per orang?

Distribusi ini membantu kita memahami seberapa besar variasi yang bisa terjadi dalam kelompok. Kita juga bisa membuat grafik dari hasil simulasi untuk melihat sebaran jumlah orang yang mengalami insomnia.

Rumus

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

dbinom(45, size = 60, prob = 0.7)  # Probabilitas 45 dari 60 mengalami insomnia

## [1] 0.08167591

Simulasi Binomial

set.seed(2)
binom_simul <- rbinom(100, size = 60, prob = 0.7)
hist(binom_simul, breaks = 10, col = "skyblue", main = "Simulasi Distribusi Binomial")

Interpretasi

Distribusi ini memperlihatkan variabilitas jumlah penderita insomnia di antara peminum kopi, dengan rata-rata mendekati $np = 60 \times 0.7 = 42$.

3.3 Distribusi Multinomial

Distribusi Multinomial digunakan ketika ada lebih dari dua kategori hasil, misalnya: “Tidak Insomnia”, “Insomnia Ringan”, dan “Insomnia Berat”.

Dengan asumsi peluang untuk masing-masing kategori, kita bisa menghitung berapa orang yang masuk ke dalam setiap kategori dari sejumlah orang tertentu. Ini sangat berguna jika kita ingin menganalisis tingkat keparahan insomnia, bukan hanya ada atau tidak.

Rumus

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

set.seed(3)
rmultinom(1, size = 60, prob = c(0.5, 0.2, 0.3))  # 60 peminum kopi

##      [,1]
## [1,]   28
## [2,]   12
## [3,]   20

Interpretasi

Dapat digunakan untuk mengestimasi distribusi derajat keparahan insomnia di antara peminum kopi.

3.4 Distribusi Poisson

Distribusi Poisson digunakan untuk menghitung jumlah kejadian insomnia dalam jangka waktu tertentu, seperti dalam satu minggu.

Jika diketahui bahwa rata-rata seseorang mengalami insomnia 4 kali seminggu, maka distribusi Poisson dapat membantu kita memahami seberapa sering insomnia terjadi di kelompok tersebut. Ini cocok untuk data dalam bentuk jumlah kejadian, bukan sekadar ada atau tidak.

Rumus

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

set.seed(4)
poisson_data <- rpois(30, lambda = 4)
table(poisson_data)

## poisson_data
##  0  1  2  3  4  5  6  7  8 10 
##  1  1  1  5  8  5  3  1  4  1

Interpretasi

Jika rata-rata kejadian insomnia adalah 4 kali seminggu untuk kelompok tertentu, maka distribusi ini cocok digunakan.

Kesimpulan Bab 3

Distribusi Bernoulli dan Binomial cocok untuk hasil biner.
Multinomial digunakan untuk klasifikasi ke dalam lebih dari dua tingkat kategori insomnia.
Distribusi Poisson berguna untuk menghitung frekuensi insomnia dalam waktu tetap.

Distribusi-distribusi ini menjadi dasar untuk model inferensial dan prediktif pada bab-bab selanjutnya.

4. Desain Sampling pada Studi Konsumsi Kopi dan Insomnia

Dalam bab ini, kita membahas secara menyeluruh bagaimana desain sampling dapat memengaruhi kualitas dan arah inferensi dari studi observasional antara konsumsi kopi dan gangguan insomnia.

4.1 Studi Prospektif vs Retrospektif

Studi Prospektif

Dalam desain ini, subjek diamati sejak awal berdasarkan eksposur (misalnya: apakah mereka mengonsumsi kopi), lalu ditelusuri ke depan untuk melihat apakah mengalami insomnia.

Contoh: - Kita mulai dari 120 orang: 60 orang peminum kopi, dan 60 bukan peminum. - Selanjutnya, kita amati selama 30 hari dan catat siapa yang mengalami insomnia.

Studi Retrospektif

Dalam pendekatan ini, kita mulai dari outcome terlebih dahulu, yaitu insomnia. Kita pilih subjek yang mengalami atau tidak mengalami insomnia, lalu ditelusuri kembali apakah mereka minum kopi.

Contoh: - Kita kumpulkan data dari 120 orang yang diketahui status insomnianya. - Kita lihat siapa di antara mereka yang punya kebiasaan minum kopi sebelumnya.

set.seed(100)
konsumsi <- c(rep(1, 60), rep(0, 60))
insomnia_prospektif <- c(rbinom(60, 1, 0.7), rbinom(60, 1, 0.4))
insomnia_retrospektif <- sample(insomnia_prospektif)
prospektif <- data.frame(Kopi = konsumsi, Insomnia = insomnia_prospektif)
retrospektif <- data.frame(Kopi = konsumsi, Insomnia = insomnia_retrospektif)

4.2 Studi Observasional

Studi observasional hanya mencatat data yang sudah terjadi, tanpa intervensi.

Karakteristik: - Tidak ada pengacakan. - Tidak mengatur siapa yang minum kopi. - Tidak memanipulasi atau mengontrol variabel.

table(prospektif$Kopi, prospektif$Insomnia)

##    
##      0  1
##   0 34 26
##   1 16 44

table(retrospektif$Kopi, retrospektif$Insomnia)

##    
##      0  1
##   0 29 31
##   1 21 39

Interpretasi: - Dari tabel kita bisa menghitung frekuensi insomnia di antara peminum dan bukan peminum kopi. - Tapi kita tidak bisa menyimpulkan kausalitas karena kemungkinan ada confounding variable.

4.3 Studi Kasus-Kontrol vs Kohort

Kasus-Kontrol

Mulai dari hasil (misalnya insomnia), kemudian identifikasi apakah individu memiliki faktor risiko (minum kopi).
Efisien untuk kasus yang jarang.
Rentan bias recall.

kasus <- subset(prospektif, Insomnia == 1)
kontrol <- subset(prospektif, Insomnia == 0)
prop_kasus <- prop.table(table(kasus$Kopi))
prop_kontrol <- prop.table(table(kontrol$Kopi))
prop_kasus

## 
##         0         1 
## 0.3714286 0.6285714

prop_kontrol

## 
##    0    1 
## 0.68 0.32

Kohort

Mulai dari eksposur (minum kopi), lalu amati outcome.
Dapat memperkirakan probabilitas dan risiko.

agregat_kohort <- aggregate(Insomnia ~ Kopi, data = prospektif, mean)
agregat_kohort

##   Kopi  Insomnia
## 1    0 0.4333333
## 2    1 0.7333333

Interpretasi: - Angka 0.7 pada kelompok peminum kopi menunjukkan 70% dari mereka mengalami insomnia. - Dapat digunakan untuk menghitung relative risk dan attributable risk di bab selanjutnya.

4.4 Tabel Perbandingan Desain Sampling

Jenis Studi	Seleksi Berdasarkan	Cocok untuk	Contoh
Prospektif	Eksposur	Estimasi Risiko	Minum kopi → Insomnia
Retrospektif	Outcome	Asosiasi Historis	Insomnia → Riwayat kopi
Kasus-Kontrol	Outcome	Rare Event Analysis	Insomnia → Lihat eksposur
Kohort	Eksposur	Temporal Causal Link	Minum kopi → Insomnia

Kesimpulan Bab 4

Desain sampling adalah fondasi penting dalam studi observasional.
Studi prospektif dan kohort lebih kuat untuk menyimpulkan hubungan sebab-akibat.
Retrospektif dan kasus-kontrol lebih ekonomis dan cocok untuk kejadian langka.
Dalam studi insomnia, kohort cocok karena insomnia bukan kejadian sangat langka dan kita ingin melihat dampak waktu setelah konsumsi kopi.

5. Tabel Kontingensi 2 × 2

Bab ini membahas distribusi peluang dan ukuran asosiasi yang dapat dihitung dari tabel kontingensi 2 × 2. Topik dibagi menjadi dua bagian utama: distribusi peluang (bersama, marginal, bersyarat) dan ukuran asosiasi (risk difference, relative risk, dan odds ratio).

5.1 Distribusi Peluang

5.1.1 Peluang Bersama

Peluang bersama adalah probabilitas dua kejadian terjadi secara bersamaan. Contoh: peluang seseorang minum kopi dan mengalami insomnia.

library(dplyr)

## Warning: package 'dplyr' was built under R version 4.3.3

## 
## Attaching package: 'dplyr'

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

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

tabel_ko <- table(prospektif$Kopi, prospektif$Insomnia)
p_bersama <- prop.table(tabel_ko)
p_bersama

##    
##             0         1
##   0 0.2833333 0.2166667
##   1 0.1333333 0.3666667

5.1.2 Peluang Marginal

Peluang marginal dihitung dari total baris atau kolom. Contoh: peluang seseorang minum kopi (tanpa memedulikan insomnia).

prop.table(tabel_ko, margin = 1)  # dalam baris (per konsumsi kopi)

##    
##             0         1
##   0 0.5666667 0.4333333
##   1 0.2666667 0.7333333

prop.table(tabel_ko, margin = 2)  # dalam kolom (per status insomnia)

##    
##             0         1
##   0 0.6800000 0.3714286
##   1 0.3200000 0.6285714

5.1.3 Peluang Bersyarat

Peluang bersyarat menunjukkan probabilitas suatu kejadian terjadi dengan syarat kejadian lain telah terjadi.

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

Contoh: peluang mengalami insomnia jika minum kopi.

p_kopi <- tabel_ko["1", ] / sum(tabel_ko["1", ])
p_tanpa <- tabel_ko["0", ] / sum(tabel_ko["0", ])
p_kopi

##         0         1 
## 0.2666667 0.7333333

p_tanpa

##         0         1 
## 0.5666667 0.4333333

Interpretasi: - Probabilitas insomnia lebih tinggi pada peminum kopi dibanding yang tidak.

5.2 Ukuran Asosiasi

5.2.1 Risk Difference (Selisih Risiko)

\[ RD = P_{eksposur} - P_{non\text{-}eksposur} \]

insomnia_kopi <- mean(prospektif$Insomnia[prospektif$Kopi == 1])
insomnia_tanpa <- mean(prospektif$Insomnia[prospektif$Kopi == 0])
rd <- insomnia_kopi - insomnia_tanpa
rd

## [1] 0.3

5.2.2 Relative Risk (Risiko Relatif)

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

rr <- insomnia_kopi / insomnia_tanpa
rr

## [1] 1.692308

5.2.3 Odds Ratio (Rasio Odds)

\[ OR = \frac{ad}{bc} \]

a <- tabel_ko["1", "1"]
b <- tabel_ko["1", "0"]
c <- tabel_ko["0", "1"]
d <- tabel_ko["0", "0"]
or <- (a * d) / (b * c)
or

## [1] 3.596154

Ringkasan Ukuran Asosiasi

data.frame(
  Ukuran = c("Risk Difference", "Relative Risk", "Odds Ratio"),
  Nilai = c(rd, rr, or)
)

##            Ukuran    Nilai
## 1 Risk Difference 0.300000
## 2   Relative Risk 1.692308
## 3      Odds Ratio 3.596154

Kesimpulan Bab 5

Peluang bersama, marginal, dan bersyarat membantu membentuk dasar pemahaman hubungan antar kategori.
Risk Difference memberi gambaran absolut peningkatan risiko akibat kopi.
Relative Risk menunjukkan proporsi peningkatan risiko.
Odds Ratio sering digunakan dalam studi kasus-kontrol.

Ukuran-ukuran ini penting dalam menilai kekuatan asosiasi antara variabel eksposur (kopi) dan outcome (insomnia).

6. Inferensi Tabel Kontingensi Dua Arah

Bab ini membahas metode inferensial yang digunakan dalam tabel kontingensi dua arah, mencakup estimasi titik dan interval, uji proporsi dan asosiasi, serta analisis residual. Studi tetap difokuskan pada hubungan antara konsumsi kopi dan insomnia.

6.1 Estimasi Titik dan Interval

6.1.1 Estimasi Titik

Estimasi titik dilakukan untuk menghitung proporsi insomnia pada masing-masing kelompok peminum dan bukan peminum kopi.

prop_kopi <- mean(prospektif$Insomnia[prospektif$Kopi == 1])
prop_nonkopi <- mean(prospektif$Insomnia[prospektif$Kopi == 0])
prop_kopi

## [1] 0.7333333

prop_nonkopi

## [1] 0.4333333

Interpretasi: - Proporsi insomnia pada peminum kopi adalah prop_kopi. - Proporsi insomnia pada bukan peminum kopi adalah prop_nonkopi.

6.1.2 Estimasi Interval

Digunakan untuk menentukan interval kepercayaan dari perbedaan proporsi dua populasi.

ci <- prop.test(x = c(sum(prospektif$Insomnia[prospektif$Kopi == 1]), 
                     sum(prospektif$Insomnia[prospektif$Kopi == 0])), 
                n = c(sum(prospektif$Kopi == 1), sum(prospektif$Kopi == 0)))
ci

## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(sum(prospektif$Insomnia[prospektif$Kopi == 1]), sum(prospektif$Insomnia[prospektif$Kopi == 0])) out of c(sum(prospektif$Kopi == 1), sum(prospektif$Kopi == 0))
## X-squared = 9.9086, df = 1, p-value = 0.001645
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  0.1152803 0.4847197
## sample estimates:
##    prop 1    prop 2 
## 0.7333333 0.4333333

Interpretasi: - Interval konfidensi yang tidak mencakup nol mengindikasikan perbedaan proporsi insomnia yang signifikan.

6.2 Uji Hipotesis

6.2.1 Uji Proporsi

Digunakan untuk menguji apakah perbedaan proporsi insomnia antara dua kelompok signifikan.

prop.test(c(sum(prospektif$Insomnia[prospektif$Kopi == 1]), 
            sum(prospektif$Insomnia[prospektif$Kopi == 0])),
          c(sum(prospektif$Kopi == 1), sum(prospektif$Kopi == 0)))

## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(sum(prospektif$Insomnia[prospektif$Kopi == 1]), sum(prospektif$Insomnia[prospektif$Kopi == 0])) out of c(sum(prospektif$Kopi == 1), sum(prospektif$Kopi == 0))
## X-squared = 9.9086, df = 1, p-value = 0.001645
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  0.1152803 0.4847197
## sample estimates:
##    prop 1    prop 2 
## 0.7333333 0.4333333

6.2.2 Uji Asosiasi (Chi-Square)

Digunakan untuk menguji apakah terdapat asosiasi antara dua variabel kategorik.

tabel_kontingensi <- table(prospektif$Kopi, prospektif$Insomnia)
chisq.test(tabel_kontingensi, correct = FALSE)

## 
##  Pearson's Chi-squared test
## 
## data:  tabel_kontingensi
## X-squared = 11.109, df = 1, p-value = 0.0008593

6.2.3 Uji Independensi (Fisher Exact)

Alternatif Chi-Square untuk sampel kecil atau sel dengan frekuensi rendah.

fisher.test(tabel_kontingensi)

## 
##  Fisher's Exact Test for Count Data
## 
## data:  tabel_kontingensi
## p-value = 0.001523
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
##  1.566297 8.352563
## sample estimates:
## odds ratio 
##   3.555572

Interpretasi: - p-value < 0.05 pada salah satu uji menandakan adanya asosiasi signifikan antara kopi dan insomnia.

6.3 Analisis Residual

6.3.1 Jenis Residual

Residual Pearson dan standardized residual digunakan untuk mengidentifikasi sel yang berkontribusi besar terhadap nilai Chi-Square.

res_pearson <- chisq.test(tabel_kontingensi)$residuals
res_std <- chisq.test(tabel_kontingensi)$stdres
res_pearson

##    
##             0         1
##   0  1.800000 -1.521278
##   1 -1.800000  1.521278

res_std

##    
##             0         1
##   0  3.332952 -3.332952
##   1 -3.332952  3.332952

6.3.2 Deteksi Outlier

Interpretasi: - Residual dengan nilai absolut > 2 menunjukkan sel yang menyimpang signifikan dari nilai harapan. - Analisis residual memberikan wawasan tentang sel mana yang paling berkontribusi terhadap asosiasi.

Kesimpulan Bab 6

Estimasi titik memberikan gambaran proporsi insomnia pada masing-masing kelompok.
Estimasi interval memberikan rentang kepercayaan terhadap perbedaan proporsi.
Uji proporsi, Chi-Square, dan Fisher memberikan validasi statistik terhadap hubungan antar variabel.
Analisis residual membantu mendeteksi penyumbang utama terhadap asosiasi.

library(dplyr)
library(ggplot2)

7. Tabel Kontingensi Tiga Arah

Bab ini membahas analisis lanjutan terhadap tiga variabel kategorik melalui tabel kontingensi tiga arah. Fokus diberikan pada eksplorasi peluang marginal, peluang bersyarat, ukuran asosiasi dalam tabel parsial, independensi bersyarat, serta pengujian statistika terkait.

7.1 Tabel Parsial dan Marginal

Kita membentuk tabel tiga arah dari data simulasi untuk tiga variabel kategorik: konsumsi kopi, insomnia, dan jenis kelamin.

set.seed(100)
kopi <- sample(c("Ya", "Tidak"), 200, replace = TRUE)
insomnia <- sample(c("Ya", "Tidak"), 200, replace = TRUE)
jenis_kelamin <- sample(c("Laki", "Perempuan"), 200, replace = TRUE)
df3 <- data.frame(Kopi = kopi, Insomnia = insomnia, Gender = jenis_kelamin)
ftab <- xtabs(~ Kopi + Insomnia + Gender, data = df3)
ftab

## , , Gender = Laki
## 
##        Insomnia
## Kopi    Tidak Ya
##   Tidak    29 26
##   Ya       21 30
## 
## , , Gender = Perempuan
## 
##        Insomnia
## Kopi    Tidak Ya
##   Tidak    25 24
##   Ya       24 21

Interpretasi: - Tabel ini menggambarkan distribusi frekuensi untuk semua kombinasi antara kopi, insomnia, dan jenis kelamin. - Dari sini kita bisa mengekstrak tabel marginal dan parsial.

Tabel marginal:

margin.table(ftab, margin = c(1, 3)) # Marginal untuk Kopi dan Gender

##        Gender
## Kopi    Laki Perempuan
##   Tidak   55        49
##   Ya      51        45

margin.table(ftab, margin = c(2, 3)) # Marginal untuk Insomnia dan Gender

##         Gender
## Insomnia Laki Perempuan
##    Tidak   50        49
##    Ya      56        45

7.2 Distribusi Peluang

Kita hitung distribusi peluang untuk ketiga variabel.

prop.table(ftab)

## , , Gender = Laki
## 
##        Insomnia
## Kopi    Tidak    Ya
##   Tidak 0.145 0.130
##   Ya    0.105 0.150
## 
## , , Gender = Perempuan
## 
##        Insomnia
## Kopi    Tidak    Ya
##   Tidak 0.125 0.120
##   Ya    0.120 0.105

Interpretasi: - Probabilitas untuk setiap sel menggambarkan peluang terjadinya kombinasi kategori secara proporsional dari total keseluruhan data.

7.3 Tabel Peluang Bersyarat

Distribusi peluang bersyarat memperlihatkan bagaimana distribusi dua variabel berubah ketika variabel ketiga dikondisikan.

prop.table(ftab, margin = 3)  # dikondisikan pada Gender

## , , Gender = Laki
## 
##        Insomnia
## Kopi        Tidak        Ya
##   Tidak 0.2735849 0.2452830
##   Ya    0.1981132 0.2830189
## 
## , , Gender = Perempuan
## 
##        Insomnia
## Kopi        Tidak        Ya
##   Tidak 0.2659574 0.2553191
##   Ya    0.2553191 0.2234043

Interpretasi: - Perbandingan pola distribusi antara kelompok Laki dan Perempuan bisa dilihat untuk mengetahui perbedaan pola asosiasi.

7.4 Ukuran Asosiasi

7.4.1 Tabel Kontingensi Parsial

Kita pecah tabel ke dalam dua strata (parsial) berdasarkan jenis kelamin dan hitung ukuran asosiasi masing-masing.

ftab_laki <- xtabs(~ Kopi + Insomnia, data = subset(df3, Gender == "Laki"))
ftab_perempuan <- xtabs(~ Kopi + Insomnia, data = subset(df3, Gender == "Perempuan"))
ftab_laki

##        Insomnia
## Kopi    Tidak Ya
##   Tidak    29 26
##   Ya       21 30

ftab_perempuan

##        Insomnia
## Kopi    Tidak Ya
##   Tidak    25 24
##   Ya       24 21

Kemudian, kita hitung Odds Ratio untuk masing-masing strata:

OR <- function(tab) {
  (tab[1,1]*tab[2,2]) / (tab[1,2]*tab[2,1])
}
OR_laki <- OR(ftab_laki)
OR_perempuan <- OR(ftab_perempuan)
OR_laki

## [1] 1.593407

OR_perempuan

## [1] 0.9114583

Interpretasi: - Nilai OR menunjukkan kekuatan asosiasi antara kopi dan insomnia pada masing-masing strata gender.

7.5 Conditional Independence

Independensi bersyarat menguji apakah dua variabel (kopi dan insomnia) independen jika dikondisikan terhadap variabel ketiga (gender).

mosaicplot(ftab, main = "Mosaic Plot Tabel Kontingensi Tiga Arah")

7.6 Marginal Y dan X

Kita ambil distribusi marginal dari insomnia (Y) dan kopi (X).

marginal_Y <- margin.table(ftab, margin = 2)
marginal_X <- margin.table(ftab, margin = 1)
marginal_Y

## Insomnia
## Tidak    Ya 
##    99   101

marginal_X

## Kopi
## Tidak    Ya 
##   104    96

Interpretasi: - Marginal ini digunakan untuk melihat distribusi total insomnia dan kopi terlepas dari gender.

7.7 Inferensi Tabel Kontingensi Tiga Arah

7.7.1 Independensi Bersyarat

Untuk menyelidiki independensi bersyarat, kita lakukan uji Chi-Square terpisah untuk masing-masing strata gender.

chisq.test(ftab_laki)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  ftab_laki
## X-squared = 0.99118, df = 1, p-value = 0.3195

chisq.test(ftab_perempuan)

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  ftab_perempuan
## X-squared = 0.00030934, df = 1, p-value = 0.986

7.7.2 Pengujian Statistik untuk Independensi Bersyarat

Uji ini akan menunjukkan apakah terdapat asosiasi signifikan antara kopi dan insomnia setelah dikondisikan pada gender.

7.7.3 Odds Ratio Bersama

Gabungkan odds ratio dari strata laki-laki dan perempuan untuk menghasilkan ukuran asosiasi umum.

OR_common <- (OR_laki + OR_perempuan) / 2
OR_common

## [1] 1.252432

Kesimpulan Bab 7

Tabel kontingensi tiga arah memperluas analisis hubungan variabel kategorik.
Peluang marginal dan bersyarat membantu memahami struktur data lebih dalam.
Pengujian OR dan independensi bersyarat memungkinkan deteksi efek confounding.
Penting memeriksa hasil pada tiap strata agar tidak keliru dalam menarik kesimpulan global.

8 Generalized Linear Model (GLM)

Generalized Linear Model (GLM) adalah perluasan dari regresi linear klasik yang memungkinkan pemodelan data dengan variabel respons yang tidak berdistribusi normal atau hubungan non-linear antara prediktor dan respons. GLM memiliki tiga komponen utama:

Distribusi exponential family untuk variabel respons.
Fungsi link yang menghubungkan ekspektasi respons dengan kombinasi linear prediktor.
Prediktor linear: $\eta = X \beta$

8.1 Exponential Family

Distribusi dalam exponential family memiliki bentuk:

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

Contoh distribusi exponential family: - Normal - Binomial - Poisson - Gamma

Contoh: Distribusi Binomial

Fungsi probabilitas binomial:

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

Dalam bentuk exponential family:

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

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

8.2 Model Regresi Logistik

Regresi logistik mirip regresi linear, tetapi membatasi output ke nilai biner (0 atau 1) menggunakan fungsi sigmoid. Output berada pada rentang 0 hingga 1, membentuk kurva S. Model ini mengklasifikasikan data ke dalam kategori diskrit berdasarkan variabel independen, sering digunakan untuk prediksi probabilitas.

Contoh Penerapan: - Memprediksi risiko insomnia berdasarkan konsumsi kopi, usia, dan gender. - Mengklasifikasi email spam. - Memperkirakan probabilitas diterima di universitas.

Keunggulan Regresi Logistik: 1. Implementasi mudah: Komputasi rendah, mudah ditafsirkan. 2. Data linear terpisah: Efektif untuk klasifikasi dua kelas yang dapat dipisahkan garis lurus. 3. Wawasan prediktor: Menunjukkan relevansi dan arah hubungan prediktor-respons.

Fungsi Sigmoid:

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

Prediksi > 0.5 diklasifikasikan sebagai 1 (positif), dan < 0.5 sebagai 0 (negatif).

Simulasi Regresi Logistik

set.seed(42)
n <- 100
coffee <- seq(0, 5, length.out = n)  # Konsumsi kopi (gelas/hari)
log_odds <- -1.2 + 0.8 * coffee
prob <- 1 / (1 + exp(-log_odds))
insomnia <- rbinom(n, 1, prob)
data <- data.frame(coffee = coffee, insomnia = insomnia, prob = prob)

# Visualisasi
plot(coffee, insomnia, pch = 16, col = "gray60",
     xlab = "Konsumsi Kopi (gelas/hari)", ylab = "Insomnia (0/1) / Probabilitas",
     main = "Regresi Logistik: Konsumsi Kopi vs Insomnia")
lines(coffee, prob, col = "blue", lwd = 2)
abline(h = 0.5, col = "red", lty = 2)
legend("topleft",
       legend = c("Data (0/1)", "Kurva Sigmoid", "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")

Visualisasi

plot(coffee, insomnia, pch = 16, col = “gray60”, xlab = “Konsumsi Kopi (gelas/hari)”, ylab = “Insomnia (0/1) / Probabilitas”, main = “Regresi Logistik: Konsumsi Kopi vs Insomnia”) lines(coffee, prob, col = “blue”, lwd = 2) abline(h = 0.5, col = “red”, lty = 2) legend(“topleft”, legend = c(“Data (0/1)”, “Kurva Sigmoid”, “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”)

Fungsi Link (Logit): \[ g(\mu) = \log \left( \frac{\mu}{1 - \mu} \right) \] Model: \[ \log \left( \frac{\mu}{1 - \mu} \right) = X \beta \] Inverse Link: \[ \mu = \frac{\exp(X \beta)}{1 + \exp(X \beta)} \] Estimasi Parameter: Log-likelihood regresi logistik: \[ l(\beta) = \sum_{i=1}^n \left[ y_i \log(\pi_i) + (1 - y_i) \log(1 - \pi_i) \right] \] dengan: \[ \pi_i = \frac{\exp(x_i^\top \beta)}{1 + \exp(x_i^\top \beta)} \]

set.seed(123)
n <- 200
coffee <- runif(n, 0, 5)  # Konsumsi kopi (gelas/hari)
age <- rnorm(n, 30, 5)    # Usia
gender <- sample(c("Male", "Female"), n, TRUE)  # Gender
log_odds <- -1.5 + 0.7 * coffee + 0.05 * age + ifelse(gender == "Male", 0.3, 0)
prob <- 1 / (1 + exp(-log_odds))
insomnia <- rbinom(n, 1, prob)
data <- data.frame(insomnia, coffee, age, gender)

# Estimasi model
model <- glm(insomnia ~ coffee + age + gender, data = data, family = binomial)
summary(model)

## 
## Call:
## glm(formula = insomnia ~ coffee + age + gender, family = binomial, 
##     data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -1.79328    1.36804  -1.311    0.190    
## coffee       0.95712    0.21155   4.524 6.06e-06 ***
## age          0.04492    0.04283   1.049    0.294    
## genderMale   0.62722    0.44251   1.417    0.156    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 169.08  on 199  degrees of freedom
## Residual deviance: 138.10  on 196  degrees of freedom
## AIC: 146.1
## 
## Number of Fisher Scoring iterations: 6

# Visualisasi
data$pred <- predict(model, type = "response")
library(ggplot2)
ggplot(data, aes(x = coffee, y = insomnia, color = gender)) +
  geom_point(alpha = 0.5) +
  geom_line(aes(y = pred), linewidth = 1.5) +
  facet_wrap(~gender) +
  labs(title = "Konsumsi Kopi vs Probabilitas Insomnia",
       x = "Konsumsi Kopi (gelas/hari)",
       y = "Insomnia (0/1) / Probabilitas",
       color = "Gender") +
  theme_minimal()

##8.3 Model Regresi Poisson Regresi Poisson digunakan untuk data cacah (bilangan bulat non-negatif) dengan distribusi Poisson.

Probabilitas Poisson: \[ P(Y = y) = \frac{e^{-\lambda} \lambda^y}{y!} \] Dalam bentuk exponential family: \[ f(y; \theta) = \exp \left( y \log(\lambda) - \lambda - \log(y!) \right) \] dengan:

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

Fungsi Link: \[ g(\mu) = \log(\mu) \] Model: \[ \log(\mu_i) = x_i^\top \beta \] Inverse Link: \[ \mu_i = \exp(x_i^\top \beta) \] Estimasi Parameter: Log-likelihood regresi Poisson: \[ l(\beta) = \sum_{i=1}^n \left[ y_i x_i^\top \beta - \exp(x_i^\top \beta) - \log(y_i!) \right] \]

#Contoh Kasus:

set.seed(42)
n <- 200
coffee <- runif(n, 0, 5)
age <- rnorm(n, 30, 5)
lambda <- exp(0.5 + 0.4 * coffee + 0.02 * age)
episodes <- rpois(n, lambda)
data <- data.frame(episodes, coffee, age)

# Estimasi model
poisson_model <- glm(episodes ~ coffee + age, data = data, family = poisson)
summary(poisson_model)

## 
## Call:
## glm(formula = episodes ~ coffee + age, family = poisson, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) 0.332179   0.154223   2.154   0.0312 *  
## coffee      0.384197   0.017709  21.695  < 2e-16 ***
## age         0.026734   0.004624   5.782 7.39e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 765.15  on 199  degrees of freedom
## Residual deviance: 202.62  on 197  degrees of freedom
## AIC: 992.1
## 
## Number of Fisher Scoring iterations: 4

# Visualisasi
data$pred <- predict(poisson_model, type = "response")
ggplot(data, aes(x = coffee, y = episodes)) +
  geom_point(color = "darkgray", alpha = 0.6) +
  geom_line(aes(y = pred), color = "blue", linewidth = 1.5) +
  labs(title = "Konsumsi Kopi vs Jumlah Episode Insomnia",
       x = "Konsumsi Kopi (gelas/hari)",
       y = "Jumlah Episode Insomnia") +
  theme_minimal()

9 Inferensi GLM

Inferensi statistik dalam Generalized Linear Model (GLM) memerlukan pemahaman tentang ekspektasi dan varians estimator untuk mendukung uji seperti Wald test, Likelihood Ratio Test (Statistik Rasio Kemungkinan), dan pembentukan interval kepercayaan.

Ekspektasi dan Varians dalam GLM

1. Ekspektasi Estimator

Ekspektasi menunjukkan apakah estimator tidak bias:

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

Estimasi Maximum Likelihood Estimation (MLE) dalam GLM bersifat asymptotically unbiased.

2. Varians Estimator

Varians mengukur presisi estimasi parameter:

\[ \operatorname{Var}(\hat{\beta}) \approx \left[ \mathbf{X}^{\top} \mathbf{W} \mathbf{X} \right]^{-1} \]

Matriks bobot $\mathbf{W}$ bergantung pada distribusi dan fungsi link. Untuk sampel besar, estimator mengikuti distribusi normal:

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

Ini digunakan untuk: - Uji Wald - Interval kepercayaan - Nilai p

Varians Tidak Konstan

Berbeda dengan regresi linear (OLS) yang mengasumsikan homoskedastisitas:

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

Dalam GLM:

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

$\phi$: Parameter dispersi
$V(\mu)$: Fungsi varians

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

Contoh Regresi Poisson: Jumlah Episode Insomnia

set.seed(123)
n <- 100
coffee <- runif(n, 0, 5)  # Konsumsi kopi (gelas/hari)
mu <- exp(0.4 + 0.7 * coffee)
y <- rpois(n, mu)
data <- data.frame(y, coffee)
model <- glm(y ~ coffee, family = poisson)
summary(model)

## 
## Call:
## glm(formula = y ~ coffee, family = poisson)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.32338    0.09566   3.381 0.000723 ***
## coffee       0.72957    0.02461  29.646  < 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: 1212.887  on 99  degrees of freedom
## Residual deviance:   65.077  on 98  degrees of freedom
## AIC: 460.88
## 
## Number of Fisher Scoring iterations: 4

Ekspektasi dan Varians:

Ekspektasi: $\mathbb{E}[Y_i] = \mu_i$ Varians: $\operatorname{Var}(Y_i) = \mu_i$

Kesimpulan:

Ekspektasi menunjukkan ketakbiasan. Varians mendukung uji statistik. Distribusi asimptotik $\beta$ bergantung pada kedua konsep ini.

9 Inferensi GLM

Ekspektasi dan Varians dalam GLM

1. Ekspektasi Estimator

Ekspektasi menunjukkan apakah estimator tidak bias:

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

Estimasi Maximum Likelihood Estimation (MLE) dalam GLM bersifat asymptotically unbiased.

2. Varians Estimator

Varians mengukur presisi estimasi parameter:

\[ \operatorname{Var}(\hat{\beta}) \approx \left[ \mathbf{X}^{\top} \mathbf{W} \mathbf{X} \right]^{-1} \]

Matriks bobot $\mathbf{W}$ bergantung pada distribusi dan fungsi link. Untuk sampel besar, estimator mengikuti distribusi normal:

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

Ini digunakan untuk: - Uji Wald - Interval kepercayaan - Nilai p

Varians Tidak Konstan

Berbeda dengan regresi linear (OLS) yang mengasumsikan homoskedastisitas:

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

Dalam GLM:

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

$\phi$: Parameter dispersi
$V(\mu)$: Fungsi varians

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

Contoh Regresi Poisson: Jumlah Episode Insomnia

set.seed(123)
n <- 100
coffee <- runif(n, 0, 5)  # Konsumsi kopi (gelas/hari)
mu <- exp(0.4 + 0.7 * coffee)
y <- rpois(n, mu)
data <- data.frame(y, coffee)
model <- glm(y ~ coffee, family = poisson)
summary(model)

## 
## Call:
## glm(formula = y ~ coffee, family = poisson)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.32338    0.09566   3.381 0.000723 ***
## coffee       0.72957    0.02461  29.646  < 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: 1212.887  on 99  degrees of freedom
## Residual deviance:   65.077  on 98  degrees of freedom
## AIC: 460.88
## 
## Number of Fisher Scoring iterations: 4

Ekspektasi dan Varians:

Ekspektasi: $\mathbb{E}[Y_i] = \mu_i$ Varians: $\operatorname{Var}(Y_i) = \mu_i$

Kesimpulan:

Ekspektasi menunjukkan ketakbiasan. Varians mendukung uji statistik. Distribusi asimptotik $\beta$ bergantung pada kedua konsep ini.

9.1 Mencari Ekspektasi dan Varians dalam GLM

Ekspektasi Untuk distribusi exponential family: \[\log f(y; \theta) = y \theta - b(\theta) + c(y)\] Score function: \[U(\theta) = \frac{\partial \ell}{\partial \theta} = y - b'(\theta)\] Ekspektasi score: \[E[U(\theta)] = \mu - b'(\theta) = 0\] Maka: \[\mu = b'(\theta)\] Varians Turunan kedua: \[\frac{\partial^2 \ell}{\partial \theta^2} = -b''(\theta)\] Sehingga: \[\operatorname{Var}(Y) = b''(\theta) = \phi V(\mu)\]

#9.2 Metode Penaksiran Parameter Maximum Likelihood Estimation (MLE)

Memaksimalkan likelihood atau log-likelihood. Syarat: Turunan pertama = 0, turunan kedua < 0.

Karena GLM tidak eksplisit, digunakan metode numerik: Newton-Raphson

Menggunakan gradien (score vector) dan matriks Hessian. Iterasi:

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

Mengganti Hessian dengan matriks informasi Fisher.

Iteratively Reweighted Least Squares (IRLS)

Mirip Fisher Scoring, serupa dengan Least Squares.

Implementasi Newton-Raphson Score function: \[U_j(\beta) = \frac{\partial \log L(\beta)}{\partial \beta_j}\] Hessian: \[H_{jk}(\beta) = \frac{\partial^2 \log L(\beta)}{\partial \beta_j \partial \beta_k}\] Iterasi: \[\hat{\beta} \approx \beta^{(t)} - H^{-1}(\beta^{(t)}) U(\beta^{(t)})\]

#9.3 Diagnostik Model GLM Diagnostik mengevaluasi kecocokan model:

Uji formal Plot prediksi vs observasi

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

Devians kecil menunjukkan model cocok.

Statistik Chi-Kuadrat Pearson: \[X^2 = \sum \frac{(y_i - \hat{\mu}_i)^2}{\hat{\mu}_i}\]

Signifikan menunjukkan model lebih baik dari model null.

Analisis Residual:

Residual = observasi - prediksi. Digunakan untuk mendeteksi penyimpangan.

9.4 Metode Estimasi dan Inferensi Regresi Logistik Regresi logistik memodelkan probabilitas respons biner:

\[\pi(x) = \frac{\exp(\beta_0 + \beta_1 x)}{1 + \exp(\beta_0 + \beta_1 x)}\] Log-likelihood: \[\ell(\beta) = \sum_{i=1}^n \left[ y_i \log(\pi_i) + (1 - y_i) \log(1 - \pi_i) \right]\] Estimasi dengan Newton-Raphson Score Function: \[\mathbf{U}(\beta) = \frac{\partial \ell(\beta)}{\partial \beta} = \mathbf{X}^{\top}(\mathbf{y} - \pi)\] Hessian Matrix: \[\mathbf{H}(\beta) = -\mathbf{X}^{\top} \mathbf{W} \mathbf{X}, \quad \text{dengan } \mathbf{W} = \operatorname{diag}(\pi_i (1 - \pi_i))\] Iterasi: \[\beta^{(t+1)} = \beta^{(t)} + \left( \mathbf{X}^{\top} \mathbf{W}^{(t)} \mathbf{X} \right)^{-1} \mathbf{X}^{\top} (\mathbf{y} - \pi^{(t)})\]

#Contoh: Konsumsi Kopi, Usia, dan Gender terhadap Insomnia

set.seed(123)
n <- 100
coffee <- runif(n, 0, 5)  # Konsumsi kopi (gelas/hari)
age <- rnorm(n, 30, 5)    # Usia
gender <- sample(c("Male", "Female"), n, TRUE)
log_odds <- -1.2 + 0.6 * coffee + 0.04 * age + ifelse(gender == "Male", 0.2, 0)
p <- 1 / (1 + exp(-log_odds))
y <- rbinom(n, 1, p)
X <- cbind(1, coffee, age, as.numeric(gender == "Male"))
data <- data.frame(y, coffee, age, gender)

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

##               [,1]
##        -2.05559048
## coffee  0.52817650
## age     0.06181858
##         0.99676787

# Menggunakan glm
model <- glm(y ~ coffee + age + gender, family = binomial, data = data)
summary(model)

## 
## Call:
## glm(formula = y ~ coffee + age + gender, family = binomial, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)  
## (Intercept) -2.05559    1.85602  -1.108   0.2681  
## coffee       0.52818    0.20779   2.542   0.0110 *
## age          0.06182    0.05833   1.060   0.2892  
## genderMale   0.99677    0.55065   1.810   0.0703 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 97.245  on 99  degrees of freedom
## Residual deviance: 86.864  on 96  degrees of freedom
## AIC: 94.864
## 
## Number of Fisher Scoring iterations: 5

9.4.1 Menghitung Statistik Rasio Kemungkinan secara Manual

Rumus: $\text{Statistik Rasio Kemungkinan} = -2 (\ell_0 - \ell_1)$$\$\ell\_0$: Log-likelihood model null: $\ell_1$

Log-likelihood: \[\ell = \sum_{i=1}^n \left[ y_i \log(p_i) + (1 - y_i) \log(1 - p_i) \right]\]

Contoh Kasus Mini: Konsumsi Kopi

data <- data.frame(
  y = c(1, 0, 1, 0),
  coffee = c(2, 1, 3, 2)
)

# Model Null
p_null <- mean(data$y)
loglik_null <- sum(data$y * log(p_null) + (1 - data$y) * log(1 - p_null))

# Model Full (misal beta0 = -1, beta1 = 1)
eta <- -1 + data$coffee
p_hat <- 1 / (1 + exp(-eta))
loglik_full <- sum(data$y * log(p_hat) + (1 - data$y) * log(1 - p_hat))

# Statistik Rasio Kemungkinan
lrt_stat <- -2 * (loglik_null - loglik_full)
p_value <- pchisq(lrt_stat, df = 1, lower.tail = FALSE)

Interpretasi:

Statistik Rasio Kemungkinan besar menunjukkan model penuh lebih baik. Nilai p dibandingkan dengan distribusi $\chi^2$ (df = jumlah prediktor).

##9.5 Detail Metode Estimasi dan Inferensi Regresi Poisson

Regresi Poisson memodelkan data cacah: \[P(Y_i = y_i) = \frac{e^{-\lambda_i} \lambda_i^{y_i}}{y_i!}\]

Model: \[\log(\lambda_i) = \mathbf{x}_i^{\top} \beta\] Log-likelihood: \[\ell(\beta) = \sum_{i=1}^n \left[ y_i \log(\lambda_i) - \lambda_i - \log(y_i!) \right]\] dengan: \[\lambda_i = \exp(\mathbf{x}_i^{\top} \beta)\] IRLS Manual untuk Regresi Poisson Iterasi: \[\beta^{(t+1)} = \left( \mathbf{X}^{\top} \mathbf{W}^{(t)} \mathbf{X} \right)^{-1} \mathbf{X}^{\top} \mathbf{W}^{(t)} \mathbf{z}^{(t)}\]

Contoh: Jumlah Episode Insomnia

set.seed(123)
n <- 100
coffee <- runif(n, 0, 5)
age <- rnorm(n, 30, 5)
lambda <- exp(0.4 + 0.5 * coffee + 0.03 * age)
y <- rpois(n, lambda)
X <- cbind(1, coffee, age)
data <- data.frame(y, coffee, age)

# IRLS Manual
beta <- c(0, 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) break
  beta <- beta_new
}
beta

##              [,1]
##        0.55398020
## coffee 0.48856003
## age    0.02587977

# Menggunakan glm
model <- glm(y ~ coffee + age, family = poisson, data = data)
summary(model)

## 
## Call:
## glm(formula = y ~ coffee + age, family = poisson, data = data)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) 0.553980   0.168915    3.28  0.00104 ** 
## coffee      0.488560   0.019810   24.66  < 2e-16 ***
## age         0.025880   0.005115    5.06 4.19e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 827.40  on 99  degrees of freedom
## Residual deviance: 108.83  on 97  degrees of freedom
## AIC: 547.07
## 
## Number of Fisher Scoring iterations: 4

Kesimpulan:

Estimasi MLE untuk regresi Poisson menggunakan IRLS. Uji Wald dan Statistik Rasio Kemungkinan mendukung pengujian hipotesis. AIC dan BIC mengevaluasi kecocokan model.

10.1 Simulasi Data

Kita akan membuat variabel baru: jumlah cangkir kopi yang dikonsumsi per hari (ordinal: 0–4 cangkir) dan status insomnia (biner).

set.seed(123)
n <- 300
kopi_cup <- sample(0:4, n, replace = TRUE, prob = c(0.15, 0.2, 0.3, 0.2, 0.15))
insomnia <- rbinom(n, 1, prob = plogis(-1 + 0.5 * kopi_cup))
data10 <- data.frame(KopiCup = kopi_cup, Insomnia = insomnia)

10.2 Eksplorasi Data

table(data10$KopiCup)

## 
##  0  1  2  3  4 
## 45 67 88 61 39

table(data10$Insomnia)

## 
##   0   1 
## 158 142

# Proporsi insomnia per jumlah kopi
aggregate(Insomnia ~ KopiCup, data = data10, mean)

##   KopiCup  Insomnia
## 1       0 0.1777778
## 2       1 0.3731343
## 3       2 0.5113636
## 4       3 0.5901639
## 5       4 0.7179487

# Visualisasi
library(scales)

## Warning: package 'scales' was built under R version 4.3.3

ggplot(data10, aes(x = as.factor(KopiCup), fill = as.factor(Insomnia))) +
  geom_bar(position = "fill") +
  labs(title = "Proporsi Insomnia per Jumlah Kopi", x = "Cangkir Kopi per Hari", y = "Proporsi", fill = "Insomnia") +
  scale_y_continuous(labels = percent) +
  theme_minimal()

10.3 Perlakuan Variabel Ordinal

10.3.1 Sebagai Nominal

Dalam pendekatan ini, jumlah cangkir dianggap sebagai kategori yang tidak berurutan.

data10$KopiFaktor <- factor(data10$KopiCup)
model_nominal <- glm(Insomnia ~ KopiFaktor, data = data10, family = binomial(link = "logit"))
summary(model_nominal)

## 
## Call:
## glm(formula = Insomnia ~ KopiFaktor, family = binomial(link = "logit"), 
##     data = data10)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.5315     0.3899  -3.928 8.57e-05 ***
## KopiFaktor1   1.0127     0.4646   2.180 0.029274 *  
## KopiFaktor2   1.5769     0.4444   3.548 0.000388 ***
## KopiFaktor3   1.8961     0.4688   4.044 5.25e-05 ***
## KopiFaktor4   2.4658     0.5279   4.671 2.99e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 415.03  on 299  degrees of freedom
## Residual deviance: 381.56  on 295  degrees of freedom
## AIC: 391.56
## 
## Number of Fisher Scoring iterations: 4

Interpretasi: - Masing-masing level kopi dibandingkan terhadap referensi (misalnya 0 cangkir). - Tidak mempertimbangkan urutan/tingkatan. - Cocok jika efek tiap level tidak linier atau tidak berurutan.

10.3.2 Sebagai Rasio

Variabel ordinal diperlakukan sebagai numerik (linear):

model_rasio <- glm(Insomnia ~ KopiCup, data = data10, family = binomial(link = "logit"))
summary(model_rasio)

## 
## Call:
## glm(formula = Insomnia ~ KopiCup, family = binomial(link = "logit"), 
##     data = data10)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -1.1966     0.2396  -4.995 5.88e-07 ***
## KopiCup       0.5561     0.1043   5.331 9.77e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 415.03  on 299  degrees of freedom
## Residual deviance: 383.19  on 298  degrees of freedom
## AIC: 387.19
## 
## Number of Fisher Scoring iterations: 4

Interpretasi: - Setiap peningkatan 1 cangkir diasumsikan menaikkan log-odds insomnia secara proporsional. - Lebih hemat parameter dan memiliki interpretasi langsung. - Cocok jika hubungan terlihat linier dari eksplorasi sebelumnya.

Perbandingan Model

AIC(model_nominal, model_rasio)

##               df      AIC
## model_nominal  5 391.5600
## model_rasio    2 387.1879

Interpretasi: - Model dengan AIC lebih kecil dianggap lebih baik. - Membandingkan model nominal dan rasio untuk melihat mana yang lebih cocok.

Visualisasi Probabilitas Prediksi

data10$pred_nom <- predict(model_nominal, type = "response")
data10$pred_rasio <- predict(model_rasio, type = "response")

library(reshape2)

## Warning: package 'reshape2' was built under R version 4.3.3

data_plot <- melt(data10[, c("KopiCup", "pred_nom", "pred_rasio")], id.vars = "KopiCup")

ggplot(data_plot, aes(x = KopiCup, y = value, color = variable)) +
  stat_summary(fun = mean, geom = "line") +
  labs(title = "Prediksi Probabilitas Insomnia", y = "Probabilitas", x = "Cangkir Kopi", color = "Model") +
  theme_minimal()

Kesimpulan Bab 10

Regresi logistik fleksibel digunakan untuk prediktor kategorik dan ordinal.
Perlakuan sebagai nominal memperlakukan setiap level sebagai terpisah.
Perlakuan sebagai numerik menghasilkan model yang lebih sederhana jika asumsi linier dipenuhi.
Visualisasi proporsi awal membantu menentukan perlakuan prediktor yang paling sesuai.
Evaluasi AIC dan grafik probabilitas bantu menentukan pendekatan terbaik secara praktis.

11. Pemilihan Model Regresi Logistik dan Evaluasi

Bab ini menjelaskan proses pemilihan model regresi logistik yang optimal serta teknik evaluasinya. Topik meliputi pendekatan confirmatory dan exploratory, stepwise selection, ROC-AUC, pseudo R-squared, dan prinsip parsimony. Evaluasi model yang tepat memungkinkan pengambilan keputusan yang lebih akurat dalam konteks klasifikasi biner, seperti kasus hubungan antara konsumsi kopi, jenis kelamin, dan kemungkinan mengalami insomnia.

11.1 Membangun Model: Confirmatory dan Exploratory

Pendekatan confirmatory didasarkan pada teori yang telah ada. Kita memasukkan variabel prediktor karena didukung oleh pengetahuan sebelumnya. Exploratory sebaliknya mengandalkan data untuk menguji kombinasi variabel yang mungkin signifikan.

set.seed(42)
df <- data.frame(
  Kopi = sample(c("Ya", "Tidak"), 300, replace = TRUE),
  Insomnia = sample(c("Ya", "Tidak"), 300, replace = TRUE),
  Gender = sample(c("Laki", "Perempuan"), 300, replace = TRUE)
)

df <- df %>% mutate(
  InsomniaBin = ifelse(Insomnia == "Ya", 1, 0),
  KopiBin = ifelse(Kopi == "Ya", 1, 0),
  GenderBin = ifelse(Gender == "Laki", 1, 0)
)

Model confirmatory:

model_confirm <- glm(InsomniaBin ~ KopiBin + GenderBin, data = df, family = binomial)
summary(model_confirm)

## 
## Call:
## glm(formula = InsomniaBin ~ KopiBin + GenderBin, family = binomial, 
##     data = df)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept)  -0.2284     0.1960  -1.165    0.244
## KopiBin       0.3189     0.2323   1.373    0.170
## GenderBin     0.2073     0.2321   0.893    0.372
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 415.83  on 299  degrees of freedom
## Residual deviance: 413.08  on 297  degrees of freedom
## AIC: 419.08
## 
## Number of Fisher Scoring iterations: 3

Model exploratory:

model_exp <- glm(InsomniaBin ~ KopiBin * GenderBin, data = df, family = binomial)
summary(model_exp)

## 
## Call:
## glm(formula = InsomniaBin ~ KopiBin * GenderBin, family = binomial, 
##     data = df)
## 
## Coefficients:
##                   Estimate Std. Error z value Pr(>|z|)  
## (Intercept)        -0.3747     0.2261  -1.657   0.0975 .
## KopiBin             0.6296     0.3292   1.912   0.0558 .
## GenderBin           0.5082     0.3236   1.571   0.1163  
## KopiBin:GenderBin  -0.6259     0.4659  -1.344   0.1791  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 415.83  on 299  degrees of freedom
## Residual deviance: 411.27  on 296  degrees of freedom
## AIC: 419.27
## 
## Number of Fisher Scoring iterations: 4

Perbandingan deviance mengukur goodness of fit dari masing-masing model:

deviance(model_confirm); deviance(model_exp)

## [1] 413.0846

## [1] 411.2742

11.2 Stepwise Selection

Stepwise selection adalah metode otomatis untuk memilih variabel terbaik dalam model menggunakan kriteria statistik. Kita mulai dari model kosong dan memperluasnya.

model_null <- glm(InsomniaBin ~ 1, data = df, family = binomial)
model_full <- glm(InsomniaBin ~ KopiBin * GenderBin, data = df, family = binomial)
model_step <- step(model_null, scope = list(lower = model_null, upper = model_full), direction = "both")

## Start:  AIC=417.83
## InsomniaBin ~ 1
## 
##             Df Deviance    AIC
## <none>           415.83 417.83
## + KopiBin    1   413.88 417.88
## + GenderBin  1   414.98 418.98

summary(model_step)

## 
## Call:
## glm(formula = InsomniaBin ~ 1, family = binomial, data = df)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept)  0.02667    0.11548   0.231    0.817
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 415.83  on 299  degrees of freedom
## Residual deviance: 415.83  on 299  degrees of freedom
## AIC: 417.83
## 
## Number of Fisher Scoring iterations: 3

Hasil stepwise selection menunjukkan bahwa interaksi antar variabel dapat signifikan jika meningkatkan kualitas model secara statistik.

11.3 Evaluasi Model: ROC dan AUC

Kurva ROC (Receiver Operating Characteristic) adalah cara visual untuk mengevaluasi performa model klasifikasi dengan melihat trade-off antara sensitivitas dan 1-spesifisitas.

model <- glm(InsomniaBin ~ KopiBin * GenderBin, data = df, family = binomial)
summary(model)

## 
## Call:
## glm(formula = InsomniaBin ~ KopiBin * GenderBin, family = binomial, 
##     data = df)
## 
## Coefficients:
##                   Estimate Std. Error z value Pr(>|z|)  
## (Intercept)        -0.3747     0.2261  -1.657   0.0975 .
## KopiBin             0.6296     0.3292   1.912   0.0558 .
## GenderBin           0.5082     0.3236   1.571   0.1163  
## KopiBin:GenderBin  -0.6259     0.4659  -1.344   0.1791  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 415.83  on 299  degrees of freedom
## Residual deviance: 411.27  on 296  degrees of freedom
## AIC: 419.27
## 
## Number of Fisher Scoring iterations: 4

df$pred <- predict(model, type = "response")
summary(df$pred)  # Harus bervariasi (tidak semua nilai sama)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  0.4074  0.4074  0.5333  0.5067  0.5342  0.5634

summary(df$pred)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  0.4074  0.4074  0.5333  0.5067  0.5342  0.5634

library(pROC)

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

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

## 
## Attaching package: 'pROC'

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

df$pred <- predict(model_step, type = "response")
df$InsomniaFactor <- factor(df$InsomniaBin, levels = c(0, 1), labels = c("No", "Yes"))
roc_obj <- roc(response = df$InsomniaFactor, predictor = df$pred, levels = c("No", "Yes"), direction = ">")

plot(roc_obj, col = "darkblue", lwd = 2, main = "Kurva ROC")
abline(a = 0, b = 1, col = "magenta", lty = 2)

auc(roc_obj)

## Area under the curve: 0.5

AUC (Area Under Curve) yang mendekati 1 menunjukkan model yang sangat baik.

11.4 Pseudo R-Squared

Karena regresi logistik tidak memiliki R² klasik seperti regresi linear, kita gunakan pseudo-R² seperti McFadden:

null_ll <- logLik(model_null)
step_ll <- logLik(model_step)
1 - (as.numeric(step_ll) / as.numeric(null_ll))

## [1] 0

Alternatif pseudo-R²:

1 - (model_step$deviance / model_step$null.deviance)

## [1] 0

11.5 Tabel Klasifikasi

Tabel klasifikasi digunakan untuk membandingkan nilai prediksi terhadap nilai aktual.

df$pred_class <- ifelse(df$pred > 0.5, 1, 0)
table(Actual = df$InsomniaBin, Predicted = df$pred_class)

##       Predicted
## Actual   1
##      0 148
##      1 152

11.6 Perbandingan Model

Kita membandingkan model berdasarkan nilai AIC (Akaike Information Criterion), nilai yang lebih kecil menunjukkan model yang lebih baik:

AIC(model_confirm, model_step, model_exp)

##               df      AIC
## model_confirm  3 419.0846
## model_step     1 417.8350
## model_exp      4 419.2742

11.7 Likelihood Ratio Test

Uji ini membandingkan dua model dengan nested structure:

anova(model_null, model_step, test = "Chisq")

## Analysis of Deviance Table
## 
## Model 1: InsomniaBin ~ 1
## Model 2: InsomniaBin ~ 1
##   Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1       299     415.83                     
## 2       299     415.83  0        0

Jika p-value < 0.05, maka model yang lebih kompleks secara signifikan lebih baik.

11.8 Prinsip Parsimony

Model terbaik tidak selalu yang paling kompleks. Parsimony menekankan efisiensi model dalam menjelaskan data dengan sesedikit mungkin parameter.

length(coef(model_confirm)); length(coef(model_exp))

## [1] 3

## [1] 4

11.9 Evaluasi Klasifikasi dan Akurasi

Menggunakan confusion matrix dari paket caret:

df$pred_factor <- as.factor(df$pred_class)
df$actual_factor <- as.factor(df$InsomniaBin)
table(Predicted = df$pred_factor, Actual = df$actual_factor)

##          Actual
## Predicted   0   1
##         1 148 152

11.9.1 Sensitivitas dan Spesifisitas

Sensitivitas adalah proporsi benar positif yang terdeteksi. Spesifisitas adalah proporsi benar negatif.

tab <- table(df$actual_factor, df$pred_factor)
 sens <- tab["1", "1"] / sum(tab["1", ])
  spec <- tab["1", "1"] / sum(tab["1", ])
  cat("Sensitivitas:", round(sens, 3), "\n")

## Sensitivitas: 1

  cat("Spesifisitas:", round(spec, 3), "\n")

## Spesifisitas: 1

11.10 Kurva ROC (Detail)

Kurva detail ROC memberikan visualisasi perubahan sensitivitas terhadap spesifisitas pada berbagai threshold.

roc_df <- data.frame(
  Specificity = rev(roc_obj$specificities),
  Sensitivity = rev(roc_obj$sensitivities)
)
ggplot(roc_df, aes(x = 1 - Specificity, y = Sensitivity)) +
  geom_line(color = "blue") +
  geom_abline(slope = 1, intercept = 0, linetype = "dashed", color = "gray") +
  labs(title = "Kurva ROC Detail", x = "1 - Spesifisitas", y = "Sensitivitas") +
  theme_minimal()

11.11 Precision-Recall Curve

Precision adalah proporsi prediksi positif yang benar. Recall sama dengan sensitivitas.

df$actual_factor <- factor(df$InsomniaBin, levels = c("0", "1"))
df$pred_factor <- factor(ifelse(df$pred >= 0.5, "1", "0"), levels = c("0", "1"))

tab <- table(Actual = df$actual_factor, Predicted = df$pred_factor)

TP <- tab["1", "1"]
FP <- tab["0", "1"]
FN <- tab["1", "0"]

precision <- TP / (TP + FP)
recall <- TP / (TP + FN)

precision

## [1] 0.5066667

recall

## [1] 1

Visualisasi precision dan recall:

pr_df <- data.frame(precision = precision, recall = recall)
ggplot(pr_df, aes(x = recall, y = precision)) +
  geom_point(size = 4, color = "darkgreen") +
  geom_segment(aes(x = 0, xend = recall, y = precision, yend = precision), linetype = "dotted") +
  geom_segment(aes(x = recall, xend = recall, y = 0, yend = precision), linetype = "dotted") +
  labs(title = "Precision vs Recall", x = "Recall", y = "Precision") +
  theme_minimal()

11.12 Pseudo R-squared Tambahan

1 - (model_step$deviance / model_step$null.deviance)

## [1] 0

Validasi tambahan dengan log-likelihood ratio:

2 * (as.numeric(logLik(model_step)) - as.numeric(logLik(model_null)))

## [1] 0

Interpretasi angka pseudo-R² dan likelihood ratio memberikan pemahaman kuantitatif akan kekuatan model dalam menjelaskan variasi data.

Kesimpulan Bab 11

Model regresi logistik dapat dibangun dengan pendekatan confirmatory berbasis teori, atau exploratory berbasis data.
Evaluasi model penting untuk memastikan akurasi dan generalisasi, menggunakan metrik seperti AIC, deviance, pseudo R², ROC-AUC, precision-recall.
Prinsip parsimony mendorong penggunaan model sederhana yang efisien.
Kurva ROC dan PR membantu memahami performa model dalam klasifikasi.
Confusion matrix, sensitivitas, spesifisitas, dan likelihood test memperkaya validasi model sebelum digunakan untuk keputusan penting.

12. Apa itu Distribusi Multinomial

Bab ini membahas regresi logistik multinomial, yang digunakan ketika variabel respon memiliki lebih dari dua kategori. Model ini memperluas regresi logistik biner, memungkinkan analisis proporsi probabilitas dari beberapa kategori yang saling eksklusif. Model ini sangat bermanfaat dalam kasus klasifikasi dengan lebih dari dua kelas, contohnya seperti dalam klasifikasi tingkat keparahan insomnia berdasarkan faktor-faktor seperti konsumsi kopi dan jenis kelamin.

12.1 Studi Kasus

Untuk memahami konsep regresi logistik multinomial, kita akan menyusun sebuah studi kasus yang sederhana namun relevan. Kita tertarik pada kemungkinan seseorang mengalami salah satu dari tiga kategori insomnia: “Tidak Insomnia”, “Insomnia Ringan”, dan “Insomnia Berat”. Dua variabel prediktor utama yang digunakan adalah konsumsi kopi dan jenis kelamin.

set.seed(123)
df_multi <- data.frame(
  Kopi = sample(c("Ya", "Tidak"), 300, replace = TRUE),
  Gender = sample(c("Laki", "Perempuan"), 300, replace = TRUE)
)
df_multi$InsomniaCat <- sample(c("Tidak", "Ringan", "Berat"), 300, replace = TRUE)
df_multi$InsomniaCat <- factor(df_multi$InsomniaCat, levels = c("Tidak", "Ringan", "Berat"))

12.2 Multinomial Logistic Regression

Multinomial logistic regression digunakan untuk mengestimasi peluang seseorang berada pada satu dari beberapa kategori variabel dependen, berdasarkan satu atau lebih variabel independen.

12.2.1 Baseline-category logit model

Model logit kategori baseline membandingkan log-odds dari kategori target terhadap kategori baseline, yang dalam kasus ini kita pilih sebagai “Tidak”. Artinya, semua interpretasi model akan menjelaskan bagaimana prediktor mengubah kemungkinan seseorang beralih dari kondisi “Tidak” insomnia ke “Ringan” atau “Berat”.

library(nnet)

## Warning: package 'nnet' was built under R version 4.3.3

model_mn <- multinom(InsomniaCat ~ Kopi + Gender, data = df_multi)

## # weights:  12 (6 variable)
## initial  value 329.583687 
## iter  10 value 327.431163
## iter  10 value 327.431163
## iter  10 value 327.431163
## final  value 327.431163 
## converged

summary(model_mn)

## Call:
## multinom(formula = InsomniaCat ~ Kopi + Gender, data = df_multi)
## 
## Coefficients:
##        (Intercept)      KopiYa GenderPerempuan
## Ringan  0.03106679 -0.30219177     0.256248152
## Berat   0.19403405 -0.04482146    -0.007795524
## 
## Std. Errors:
##        (Intercept)    KopiYa GenderPerempuan
## Ringan   0.2613084 0.2927838       0.2927905
## Berat    0.2513994 0.2817005       0.2811700
## 
## Residual Deviance: 654.8623 
## AIC: 666.8623

12.2.2 Estimasi Parameter

Untuk mengevaluasi apakah suatu koefisien signifikan secara statistik, kita dapat menghitung nilai z-score dan p-value. Nilai ini berguna untuk menilai apakah pengaruh prediktor terhadap outcome dalam masing-masing kategori dapat dianggap nyata secara statistik.

z <- summary(model_mn)$coefficients / summary(model_mn)$standard.errors
pvals <- (1 - pnorm(abs(z))) * 2
pvals

##        (Intercept)   KopiYa GenderPerempuan
## Ringan   0.9053630 0.302010       0.3814689
## Berat    0.4402235 0.873582       0.9778812

12.3 Contoh Kasus

Dalam kehidupan nyata, klasifikasi insomnia menjadi tiga tingkat keparahan adalah contoh nyata penerapan regresi logistik multinomial. Variabel prediktor seperti konsumsi kafein dan jenis kelamin diketahui memengaruhi kualitas tidur dan oleh karena itu menjadi input penting dalam analisis ini.

12.4 Simulasi Data

Kita dapat memvisualisasikan frekuensi kategori insomnia berdasarkan kebiasaan minum kopi. Ini membantu kita memahami apakah terjadi perbedaan distribusi kategori insomnia antar kelompok peminum kopi dan bukan peminum kopi.

table(df_multi$InsomniaCat)

## 
##  Tidak Ringan  Berat 
##     94     95    111

ggplot(df_multi, aes(x = InsomniaCat, fill = Kopi)) +
  geom_bar(position = "dodge") +
  labs(title = "Distribusi Kategori Insomnia berdasarkan Kopi")

12.5 Estimasi Model

Model yang mempertimbangkan interaksi antara jenis kelamin dan kebiasaan minum kopi dapat memberikan informasi tambahan tentang apakah efek kopi terhadap insomnia berbeda pada laki-laki dan perempuan.

model_mn2 <- multinom(InsomniaCat ~ Kopi * Gender, data = df_multi)

## # weights:  15 (8 variable)
## initial  value 329.583687 
## iter  10 value 326.615569
## final  value 326.600352 
## converged

summary(model_mn2)

## Call:
## multinom(formula = InsomniaCat ~ Kopi * Gender, data = df_multi)
## 
## Coefficients:
##        (Intercept)     KopiYa GenderPerempuan KopiYa:GenderPerempuan
## Ringan   0.2336154 -0.6903744      -0.1158313              0.7469445
## Berat    0.2744373 -0.1791281      -0.1566532              0.2749183
## 
## Std. Errors:
##        (Intercept)    KopiYa GenderPerempuan KopiYa:GenderPerempuan
## Ringan   0.3070802 0.4245741       0.4159351              0.5887155
## Berat    0.3043544 0.3953074       0.4139267              0.5666148
## 
## Residual Deviance: 653.2007 
## AIC: 669.2007

12.6 Nilai P-Value dan Interpretasi

Nilai p-value dari setiap koefisien membantu kita menentukan apakah hubungan yang diamati pada data cukup kuat untuk dianggap signifikan. Jika nilai p lebih kecil dari 0.05, maka pengaruh prediktor dianggap signifikan terhadap outcome tersebut.

library(knitr)

## Warning: package 'knitr' was built under R version 4.3.3

z2 <- summary(model_mn2)$coefficients / summary(model_mn2)$standard.errors
pvals2 <- (1 - pnorm(abs(z2))) * 2
kable(pvals2, caption = "P-Value dari Koefisien Model Multinomial")

P-Value dari Koefisien Model Multinomial
	(Intercept)	KopiYa	GenderPerempuan	KopiYa:GenderPerempuan
Ringan	0.4467983	0.1039411	0.7806408	0.2045232
Berat	0.3672147	0.6504507	0.7050916	0.6275385

12.7 Prediksi dan Validasi

Akurasi prediksi dari model multinomial dapat diuji dengan confusion matrix. Matriks ini membandingkan hasil prediksi model dengan nilai aktual dari data.

predicted_class <- predict(model_mn2)
table(Actual = df_multi$InsomniaCat, Predicted = predicted_class)

##         Predicted
## Actual   Tidak Ringan Berat
##   Tidak      0     11    83
##   Ringan     0     13    82
##   Berat      0     13    98

12.8 Penutup

Model logistik multinomial adalah alat yang sangat kuat dalam analisis data kategorik dengan banyak kelas. Model ini sangat cocok ketika data tidak dapat disederhanakan menjadi dua kategori saja. Dalam kasus kita, model ini memungkinkan analisis mendalam terhadap tingkat keparahan insomnia yang dipengaruhi oleh kebiasaan konsumsi kopi dan faktor demografis.

12.8.1 Contoh Kasus 2 di R

Untuk memperkaya analisis, kita dapat menambahkan prediktor tambahan seperti usia dan tingkat stres. Ini memungkinkan kita mengevaluasi lebih dalam interaksi kompleks antar variabel.

df_multi2 <- df_multi %>% mutate(
  Usia = sample(20:50, 300, replace = TRUE),
  Stres = sample(c("Tinggi", "Rendah"), 300, replace = TRUE)
)
model_case2 <- multinom(InsomniaCat ~ Kopi + Gender + Usia + Stres, data = df_multi2)

## # weights:  18 (10 variable)
## initial  value 329.583687 
## iter  10 value 326.625382
## final  value 326.590466 
## converged

summary(model_case2)

## Call:
## multinom(formula = InsomniaCat ~ Kopi + Gender + Usia + Stres, 
##     data = df_multi2)
## 
## Coefficients:
##        (Intercept)      KopiYa GenderPerempuan        Usia StresTinggi
## Ringan  -0.2050695 -0.32328841      0.27746444 0.004342547   0.1723433
## Berat   -0.3707481 -0.08466115      0.02937826 0.011849179   0.2956578
## 
## Std. Errors:
##        (Intercept)    KopiYa GenderPerempuan       Usia StresTinggi
## Ringan   0.6348447 0.2946411       0.2948649 0.01620706   0.2958165
## Berat    0.6159464 0.2842688       0.2839277 0.01559244   0.2851126
## 
## Residual Deviance: 653.1809 
## AIC: 673.1809

13. Regresi Logistik Ordinal pada Kasus Insomnia dan Kepuasan Tidur

Model regresi logistik ordinal digunakan ketika variabel respons memiliki urutan, tetapi jarak antar kategorinya tidak diketahui atau tidak sama. Dalam kasus ini, tingkat kepuasan terhadap kualitas tidur dapat dikelompokkan menjadi tiga kategori: “Tidak Puas”, “Cukup Puas”, dan “Sangat Puas”. Model ini membantu menjelaskan faktor-faktor yang mempengaruhi perbedaan tingkat kepuasan tidur, terutama konsumsi kopi dan jenis kelamin.

13.1 Konsep Cumulative Logit Model

Cumulative Logit Model atau model logit kumulatif digunakan ketika data respons bersifat ordinal. Model ini menghitung probabilitas kumulatif berada pada atau di bawah kategori tertentu. Secara matematis:

\[ \log\left(\frac{P(Y \leq j)}{P(Y > j)}\right) = \theta_j - X\beta \]

di mana: - $\theta_j$: intercept spesifik untuk cutpoint j - $X$: matriks prediktor - $\beta$: vektor koefisien regresi

Model ini mengasumsikan efek prediktor tetap sama di setiap cutpoint (proportional odds).

13.2 Interpretasi Koefisien

Koefisien yang diperoleh dari model ini memberikan informasi mengenai arah dan kekuatan asosiasi antara prediktor dengan kategori kepuasan tidur. Misalnya, koefisien positif untuk variabel “Kopi” berarti bahwa konsumsi kopi cenderung meningkatkan peluang individu berada dalam kategori kepuasan tidur yang lebih tinggi (misal dari “Tidak Puas” ke “Cukup Puas”).

Namun, interpretasi kuantitatifnya berada dalam konteks log-odds kumulatif. Oleh karena itu, transformasi menjadi odds ratio sering kali diperlukan untuk interpretasi yang lebih intuitif.

13.3 Contoh Data: Kepuasan Tidur

Simulasi data dilakukan untuk mencerminkan situasi nyata: pengaruh konsumsi kopi terhadap kepuasan tidur berdasarkan jenis kelamin.

set.seed(13)
df13 <- data.frame(
  Kopi = sample(c("Ya", "Tidak"), 500, replace = TRUE),
  Gender = sample(c("Laki", "Perempuan"), 500, replace = TRUE),
  Kepuasan = ordered(sample(c("Tidak Puas", "Cukup Puas", "Sangat Puas"),
                            500, replace = TRUE, prob = c(0.3, 0.5, 0.2)),
                     levels = c("Tidak Puas", "Cukup Puas", "Sangat Puas"))
)
head(df13)

##    Kopi    Gender    Kepuasan
## 1 Tidak      Laki Sangat Puas
## 2    Ya      Laki  Cukup Puas
## 3 Tidak Perempuan  Cukup Puas
## 4    Ya      Laki Sangat Puas
## 5 Tidak      Laki  Cukup Puas
## 6    Ya      Laki Sangat Puas

Distribusi awal dapat divisualisasikan:

ggplot(df13, aes(x = Kepuasan, fill = Kopi)) +
  geom_bar(position = "fill") +
  facet_wrap(~Gender) +
  labs(y = "Proporsi", title = "Distribusi Kepuasan Tidur berdasarkan Kopi dan Gender")

13.4 Estimasi Model Ordinal

Gunakan model cumulative logit dengan link logit untuk memodelkan variabel ordinal.

library(VGAM)

## Warning: package 'VGAM' was built under R version 4.3.3

## Loading required package: stats4

## Loading required package: splines

model_ord <- vglm(Kepuasan ~ Kopi + Gender, family = cumulative(parallel = TRUE, link = "logitlink"), data = df13)
summary(model_ord)

## 
## Call:
## vglm(formula = Kepuasan ~ Kopi + Gender, family = cumulative(parallel = TRUE, 
##     link = "logitlink"), data = df13)
## 
## Coefficients: 
##                 Estimate Std. Error z value Pr(>|z|)    
## (Intercept):1   -1.04127    0.16293  -6.391 1.65e-10 ***
## (Intercept):2    1.16810    0.16506   7.077 1.47e-12 ***
## KopiYa           0.02904    0.16923   0.172    0.864    
## GenderPerempuan  0.20580    0.16926   1.216    0.224    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Names of linear predictors: logitlink(P[Y<=1]), logitlink(P[Y<=2])
## 
## Residual deviance: 1034.699 on 996 degrees of freedom
## 
## Log-likelihood: -517.3495 on 996 degrees of freedom
## 
## Number of Fisher scoring iterations: 3 
## 
## No Hauck-Donner effect found in any of the estimates
## 
## 
## Exponentiated coefficients:
##          KopiYa GenderPerempuan 
##        1.029467        1.228511

13.5 Nilai P-Value

Nilai p-value digunakan untuk mengevaluasi signifikansi statistik dari setiap variabel prediktor.

round(coef(summary(model_ord)), 4)[, 4]

##   (Intercept):1   (Intercept):2          KopiYa GenderPerempuan 
##          0.0000          0.0000          0.8637          0.2240

Jika p-value < 0.05 maka efek prediktor dianggap signifikan.

13.6 Prediksi Probabilitas

Hitung probabilitas prediksi untuk tiap kategori respon:

pred_probs <- predict(model_ord, type = "response")
head(pred_probs)

##   Tidak Puas Cukup Puas Sangat Puas
## 1  0.2609048  0.5018964   0.2371988
## 2  0.2665437  0.5014720   0.2319844
## 3  0.3024901  0.4955197   0.2019902
## 4  0.2665437  0.5014720   0.2319844
## 5  0.2609048  0.5018964   0.2371988
## 6  0.2665437  0.5014720   0.2319844

Gabungkan ke data asli:

df13 <- cbind(df13, as.data.frame(pred_probs))

Visualisasi distribusi probabilitas:

library(tidyr)

## 
## Attaching package: 'tidyr'

## The following object is masked from 'package:reshape2':
## 
##     smiths

df_long <- df13 %>%
  pivot_longer(cols = c("Tidak Puas", "Cukup Puas", "Sangat Puas"), names_to = "Kategori", values_to = "Probabilitas")

ggplot(df_long, aes(x = Kategori, y = Probabilitas, fill = Kategori)) +
  geom_boxplot() +
  ggtitle("Distribusi Probabilitas Prediksi per Kategori")

13.7 Goodness-of-Fit dan Proportional Odds

Evaluasi log-likelihood dan asumsi model:

logLik(model_ord)

## [1] -517.3495

Cek asumsi proporsional odds dengan membandingkan model paralel dan non-paralel:

model_nonparallel <- vglm(Kepuasan ~ Kopi + Gender, family = cumulative(parallel = FALSE), data = df13)
anova(model_ord, model_nonparallel, type = "I")

## Analysis of Deviance Table
## 
## Model 1: Kepuasan ~ Kopi + Gender
## Model 2: Kepuasan ~ Kopi + Gender
##   Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1       996     1034.7                     
## 2       994     1034.7  2 0.037656   0.9813

13.8 Alternatif Model Ordinal

Ketika asumsi proportional odds tidak terpenuhi, kita dapat menggunakan: - Partial Proportional Odds Model - Generalized Ordered Logit Model - Adjacent Category Logit Model

Namun perlu pertimbangan tambahan terhadap kompleksitas interpretasi dan struktur data.

13.9 Evaluasi Prediksi dan Confusion Matrix

Klasifikasi berdasarkan kategori dengan probabilitas maksimum:

df13$pred_class <- apply(pred_probs, 1, function(x) colnames(pred_probs)[which.max(x)])
tab_conf <- table(Predicted = df13$pred_class, Actual = df13$Kepuasan)
tab_conf

##             Actual
## Predicted    Tidak Puas Cukup Puas Sangat Puas
##   Cukup Puas        143        249         108

Hitung akurasi:

accuracy <- sum(diag(tab_conf)) / sum(tab_conf)
accuracy

## [1] 0.286

13.10 Asumsi Paralelisme dalam Regresi Logistik Ordinal

Asumsi paralelisme menyatakan bahwa semua prediktor mempengaruhi log-odds antar kategori secara konsisten di seluruh tingkatan. Jika asumsi ini dilanggar, maka hasil model menjadi bias.

Model logistik ordinal menawarkan pendekatan fleksibel dan informatif dalam menganalisis data dengan urutan kategori. Pada studi insomnia dan kepuasan tidur, kita dapat menyimpulkan bahwa kebiasaan minum kopi dan gender memiliki pengaruh terhadap persepsi kualitas tidur, terutama saat dianalisis dalam skala ordinal.

Visualisasi, uji asumsi, dan evaluasi prediksi memberikan gambaran lengkap tentang bagaimana model ini bekerja dan kapan ia sebaiknya diterapkan dibandingkan metode lain seperti regresi logistik biner.

14. Model Log-Linear Tiga Arah pada Kasus Insomnia, Kopi, dan Gender

Model log-linear digunakan untuk mengevaluasi hubungan tiga arah antara insomnia, konsumsi kopi, dan jenis kelamin. Pendekatan ini membantu memetakan sejauh mana interaksi antar variabel memengaruhi kecenderungan insomnia.

14.1 Tabel Kontingensi dan Model Loglinier

Kita simulasikan 300 data responden berdasarkan kebiasaan konsumsi kopi, kondisi insomnia, dan jenis kelamin.

set.seed(123)
df_log <- data.frame(
  Kopi = sample(c("Ya", "Tidak"), 300, replace = TRUE),
  Insomnia = sample(c("Ya", "Tidak"), 300, replace = TRUE),
  Gender = sample(c("Laki", "Perempuan"), 300, replace = TRUE)
)

Selanjutnya kita buat tabel kontingensi tiga arah dari variabel tersebut.

tab3d <- xtabs(~ Kopi + Insomnia + Gender, data = df_log)
tab3d

## , , Gender = Laki
## 
##        Insomnia
## Kopi    Tidak Ya
##   Tidak    36 41
##   Ya       37 40
## 
## , , Gender = Perempuan
## 
##        Insomnia
## Kopi    Tidak Ya
##   Tidak    42 27
##   Ya       35 42

ftable(tab3d)

##                Gender Laki Perempuan
## Kopi  Insomnia                      
## Tidak Tidak             36        42
##       Ya                41        27
## Ya    Tidak             37        35
##       Ya                40        42

Tabel ini menampilkan distribusi gabungan dari ketiga variabel, yang akan dianalisis dalam model log-linear.

14.2 Model Saturated

Model saturated mengasumsikan adanya semua efek utama dan interaksi dua serta tiga arah.

library(MASS)

## 
## Attaching package: 'MASS'

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

library(dplyr)
model_saturated <- loglm(~ Kopi * Insomnia * Gender, data = tab3d)
summary(model_saturated)

## Formula:
## ~Kopi * Insomnia * Gender
## attr(,"variables")
## list(Kopi, Insomnia, Gender)
## attr(,"factors")
##          Kopi Insomnia Gender Kopi:Insomnia Kopi:Gender Insomnia:Gender
## Kopi        1        0      0             1           1               0
## Insomnia    0        1      0             1           0               1
## Gender      0        0      1             0           1               1
##          Kopi:Insomnia:Gender
## Kopi                        1
## Insomnia                    1
## Gender                      1
## attr(,"term.labels")
## [1] "Kopi"                 "Insomnia"             "Gender"              
## [4] "Kopi:Insomnia"        "Kopi:Gender"          "Insomnia:Gender"     
## [7] "Kopi:Insomnia:Gender"
## attr(,"order")
## [1] 1 1 1 2 2 2 3
## 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

Deviance sangat kecil menunjukkan bahwa model ini sangat cocok, tetapi bisa terlalu kompleks.

14.3 Model Independent

Model independen mengabaikan semua interaksi dan hanya mempertimbangkan efek utama.

model_indep <- loglm(~ Kopi + Insomnia + Gender, data = tab3d)
summary(model_indep)

## Formula:
## ~Kopi + Insomnia + Gender
## attr(,"variables")
## list(Kopi, Insomnia, Gender)
## attr(,"factors")
##          Kopi Insomnia Gender
## Kopi        1        0      0
## Insomnia    0        1      0
## Gender      0        0      1
## attr(,"term.labels")
## [1] "Kopi"     "Insomnia" "Gender"  
## attr(,"order")
## [1] 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 4.591329  4 0.3318552
## Pearson          4.493034  4 0.3433742

Perbandingan deviance dengan saturated model menunjukkan seberapa buruk asumsi independen.

anova(model_indep, model_saturated)

## LR tests for hierarchical log-linear models
## 
## Model 1:
##  ~Kopi + Insomnia + Gender 
## Model 2:
##  ~Kopi * Insomnia * Gender 
## 
##           Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1   4.591329  4                                    
## Model 2   0.000000  0   4.591329         4        0.33186
## Saturated 0.000000  0   0.000000         0        1.00000

14.4 Odds Ratio dan Interpretasi

Untuk memahami asosiasi spesifik, kita hitung odds ratio berdasarkan strata jenis kelamin.

tab_laki <- xtabs(~ Kopi + Insomnia, data = df_log[df_log$Gender == "Laki",])
tab_perempuan <- xtabs(~ Kopi + Insomnia, data = df_log[df_log$Gender == "Perempuan",])

OR_laki <- (tab_laki[1,1] * tab_laki[2,2]) / (tab_laki[1,2] * tab_laki[2,1])
OR_perempuan <- (tab_perempuan[1,1] * tab_perempuan[2,2]) / (tab_perempuan[1,2] * tab_perempuan[2,1])
OR_laki; OR_perempuan

## [1] 0.9492419

## [1] 1.866667

Visualisasi odds ratio:

barplot(c(OR_laki, OR_perempuan), beside = TRUE, col = c("skyblue", "pink"),
        names.arg = c("Laki", "Perempuan"), main = "Odds Ratio Kopi dan Insomnia")
abline(h = 1, lty = 2)

14.5 Estimasi Parameter

Model parsial digunakan untuk melihat parameter dua arah saja.

model_partial <- loglm(~ Kopi * Insomnia + Kopi * Gender + Insomnia * Gender, data = tab3d)
summary(model_partial)

## Formula:
## ~Kopi * Insomnia + Kopi * Gender + Insomnia * Gender
## attr(,"variables")
## list(Kopi, Insomnia, Gender)
## attr(,"factors")
##          Kopi Insomnia Gender Kopi:Insomnia Kopi:Gender Insomnia:Gender
## Kopi        1        0      0             1           1               0
## Insomnia    0        1      0             1           0               1
## Gender      0        0      1             0           1               1
## attr(,"term.labels")
## [1] "Kopi"            "Insomnia"        "Gender"          "Kopi:Insomnia"  
## [5] "Kopi:Gender"     "Insomnia:Gender"
## attr(,"order")
## [1] 1 1 1 2 2 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
## 
## Statistics:
##                       X^2 df  P(> X^2)
## Likelihood Ratio 2.112031  1 0.1461453
## Pearson          2.108532  1 0.1464798

Model ini cukup sering digunakan jika interaksi tiga arah tidak signifikan.

14.6 Model Lebih Sederhana dan Perbandingan Model

Kita coba model lebih sederhana dan bandingkan dengan model sebelumnya.

model_simple <- loglm(~ Kopi + Insomnia + Gender + Kopi:Insomnia, data = tab3d)
summary(model_simple)

## Formula:
## ~Kopi + Insomnia + Gender + Kopi:Insomnia
## attr(,"variables")
## list(Kopi, Insomnia, Gender)
## attr(,"factors")
##          Kopi Insomnia Gender Kopi:Insomnia
## Kopi        1        0      0             1
## Insomnia    0        1      0             1
## Gender      0        0      1             0
## attr(,"term.labels")
## [1] "Kopi"          "Insomnia"      "Gender"        "Kopi:Insomnia"
## attr(,"order")
## [1] 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 3.256053  3 0.3538062
## Pearson          3.237196  3 0.3564799

anova(model_simple, model_saturated)

## LR tests for hierarchical log-linear models
## 
## Model 1:
##  ~Kopi + Insomnia + Gender + Kopi:Insomnia 
## Model 2:
##  ~Kopi * Insomnia * Gender 
## 
##           Deviance df Delta(Dev) Delta(df) P(> Delta(Dev)
## Model 1   3.256053  3                                    
## Model 2   0.000000  0   3.256053         3        0.35381
## Saturated 0.000000  0   0.000000         0        1.00000

aic_manual <- function(model) {
  dev <- model$deviance
  k <- model$df.residual - model$df  # jumlah parameter = total df - residual df
  return(dev + 2 * k)
}

aic_manual(model_simple)

## numeric(0)

aic_manual(model_saturated)

## numeric(0)

AIC digunakan untuk menilai efisiensi model. Model dengan AIC lebih rendah dianggap lebih baik.

14.7 Studi Kasus Lanjutan: Simulasi Alternatif

Kita perkuat analisis dengan studi alternatif dari simulasi data insomnia berdasarkan kelompok usia.

set.seed(999)
df_alt <- data.frame(
  Kopi = sample(c("Ya", "Tidak"), 300, replace = TRUE),
  Insomnia = sample(c("Ya", "Tidak"), 300, replace = TRUE),
  Umur = sample(c("<30", "30-50", ">50"), 300, replace = TRUE)
)
tab_alt <- xtabs(~ Kopi + Insomnia + Umur, data = df_alt)

model_alt <- loglm(~ Kopi * Insomnia + Kopi * Umur + Insomnia * Umur, data = tab_alt)
summary(model_alt)

## Formula:
## ~Kopi * Insomnia + Kopi * Umur + Insomnia * Umur
## attr(,"variables")
## list(Kopi, Insomnia, Umur)
## attr(,"factors")
##          Kopi Insomnia Umur Kopi:Insomnia Kopi:Umur Insomnia:Umur
## Kopi        1        0    0             1         1             0
## Insomnia    0        1    0             1         0             1
## Umur        0        0    1             0         1             1
## attr(,"term.labels")
## [1] "Kopi"          "Insomnia"      "Umur"          "Kopi:Insomnia"
## [5] "Kopi:Umur"     "Insomnia:Umur"
## attr(,"order")
## [1] 1 1 1 2 2 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
## 
## Statistics:
##                        X^2 df  P(> X^2)
## Likelihood Ratio 0.6907761  2 0.7079456
## Pearson          0.6904786  2 0.7080509

Model log-linear menunjukkan hubungan kompleks antara insomnia, kebiasaan konsumsi kopi, dan karakteristik lain seperti gender atau umur. Ini memungkinkan kita memahami struktur asosiasi yang tidak terlihat dari analisis dua arah sederhana. title: “Ebook ADK 2025” author: “Kevin Jonathan” date: “2025-06-24” output: html_document —

R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

summary(cars)

##      speed           dist       
##  Min.   : 4.0   Min.   :  2.00  
##  1st Qu.:12.0   1st Qu.: 26.00  
##  Median :15.0   Median : 36.00  
##  Mean   :15.4   Mean   : 42.98  
##  3rd Qu.:19.0   3rd Qu.: 56.00  
##  Max.   :25.0   Max.   :120.00

Including Plots

You can also embed plots, for example:

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.

Ebook Analisi Data Kategori

Kevin Jonathan

2025-06-11

1. Pendahuluan - Analisis Data Kategorik pada Kasus Insomnia

1.1 Latar Belakang Masalah

1.2 Tujuan Penelitian

1.3 Pentingnya Analisis Kategorik

Melalui pendekatan ini, kita tidak hanya dapat melihat hubungan antara dua variabel, tetapi juga mengontrol efek variabel lain dan menguji interaksi antar variabel secara menyeluruh.

1.4 Struktur Buku

1.5 Contoh Simulasi Data Awal

1.6 Visualisasi Awal Data

2. Metode dalam Analisis Data Kategori

2.1 Tabel Kontingensi dan Uji Chi-Square

Contoh Data (Hipotetik)

Rumus Chi-Square

Uji Statistik di R

Interpretasi

2.2 Regresi Logistik

Rumus Model

Odds Ratio

Interpretasi

2.3 Correspondence Analysis (CA)

2.4 Pohon Keputusan (Decision Tree)

3. Distribusi Probabilitas dalam Data Kategori

3.1 Distribusi Bernoulli

Rumus

Interpretasi

3.2 Distribusi Binomial

Rumus

Simulasi Binomial

Interpretasi

3.3 Distribusi Multinomial

Rumus

Interpretasi

3.4 Distribusi Poisson

Rumus

Interpretasi

Kesimpulan Bab 3

4. Desain Sampling pada Studi Konsumsi Kopi dan Insomnia

4.1 Studi Prospektif vs Retrospektif

Studi Prospektif

Studi Retrospektif

4.2 Studi Observasional

4.3 Studi Kasus-Kontrol vs Kohort

Kasus-Kontrol

Kohort

4.4 Tabel Perbandingan Desain Sampling

Kesimpulan Bab 4

5. Tabel Kontingensi 2 × 2

5.1 Distribusi Peluang

5.1.1 Peluang Bersama

5.1.2 Peluang Marginal

5.1.3 Peluang Bersyarat

5.2 Ukuran Asosiasi

5.2.1 Risk Difference (Selisih Risiko)

5.2.2 Relative Risk (Risiko Relatif)

5.2.3 Odds Ratio (Rasio Odds)

Ringkasan Ukuran Asosiasi

Kesimpulan Bab 5

6. Inferensi Tabel Kontingensi Dua Arah

6.1 Estimasi Titik dan Interval

6.1.1 Estimasi Titik

6.1.2 Estimasi Interval

6.2 Uji Hipotesis

6.2.1 Uji Proporsi

6.2.2 Uji Asosiasi (Chi-Square)

6.2.3 Uji Independensi (Fisher Exact)

6.3 Analisis Residual

6.3.1 Jenis Residual

6.3.2 Deteksi Outlier

Kesimpulan Bab 6

7. Tabel Kontingensi Tiga Arah

7.1 Tabel Parsial dan Marginal

7.2 Distribusi Peluang

7.3 Tabel Peluang Bersyarat

7.4 Ukuran Asosiasi

7.4.1 Tabel Kontingensi Parsial

7.5 Conditional Independence

7.6 Marginal Y dan X

7.7 Inferensi Tabel Kontingensi Tiga Arah