📋

Mata Kuliah: Analisis Data Kategori

Program Studi S1 Statistika · FMIPA Universitas Padjadjaran · Dosen: Prof. I Gede Nyoman Mindra Jaya, Ph.D

Kumpulan tugas mata kuliah Analisis Data Kategori yang mencakup konsep dasar data kategori, tabel kontingensi, ukuran asosiasi, hingga inferensi statistik. Seluruh analisis dilakukan menggunakan R/RStudio dengan pendekatan yang reprodusibel.

Tugas 1 — Konsep Dasar Analisis Data Kategori

📅 10 Maret 2026

1. Definisi Analisis Data Kategori

Analisis data kategori adalah cabang statistika yang digunakan untuk menganalisis variabel yang berbentuk kategori atau klasifikasi. Berbeda dengan data kontinu, data kategori merepresentasikan keanggotaan suatu objek pada kelompok tertentu dan tidak memiliki makna numerik langsung.

Karakteristik Variabel Kategori:

Variabel kategori tidak memiliki makna numerik langsung, melainkan menunjukkan keanggotaan suatu objek pada kelompok tertentu.
Terdapat dua jenis utama:
- Nominal — tidak memiliki urutan alami (contoh: jenis kelamin, golongan darah)
- Ordinal — memiliki urutan tertentu (contoh: tingkat pendidikan, tingkat kepuasan)

Contoh Penerapan:

Dalam penelitian kesehatan, metode ini sering digunakan untuk melihat hubungan antara kebiasaan gaya hidup (seperti merokok) dengan status penyakit (seperti kanker). Selain itu, bidang sosial menggunakannya untuk menganalisis tingkat kepuasan atau status ekonomi masyarakat.

2. Tabel Kontingensi

2.1 Definisi

Tabel kontingensi adalah tabel yang menyajikan distribusi frekuensi dari dua atau lebih variabel kategori secara bersamaan untuk melihat hubungan antar kategori tersebut.

2.2 Struktur Tabel Kontingensi 2×2

Contoh struktur hubungan antara variabel baris (X) dan kolom (Y):

	Kolom 1 (\(j=1\))	Kolom 2 (\(j=2\))	Total
Baris 1 (\(i=1\))	\(n_{11}\)	\(n_{12}\)	\(n_{1.}\)
Baris 2 (\(i=2\))	\(n_{21}\)	\(n_{22}\)	\(n_{2.}\)
Total	\(n_{.1}\)	\(n_{.2}\)	\(n\)

Keterangan Notasi:

\(n_{ij}\) : Jumlah observasi pada baris \(i\) dan kolom \(j\)
\(n_{i.}\) : Total marginal baris ke-\(i\)
\(n_{.j}\) : Total marginal kolom ke-\(j\)
\(n\) : Total keseluruhan sampel

2.3 Konsep Distribusi

Berdasarkan struktur tabel di atas, terdapat tiga konsep probabilitas utama:

Joint Distribution — probabilitas elemen baris \(i\) dan kolom \(j\) terjadi bersamaan: \[P_{ij} = \frac{n_{ij}}{n}\]

Marginal Distribution — probabilitas satu variabel tanpa memperhatikan variabel lain: \[P_{i.} = \sum_{j} P_{ij} \qquad P_{.j} = \sum_{i} P_{ij}\]

Conditional Probability — probabilitas suatu kejadian dengan syarat kategori tertentu sudah diketahui: \[P(j \mid i) = \frac{P_{ij}}{P_{i.}}\]

3. Ukuran Asosiasi

Ukuran asosiasi digunakan untuk mengukur kekuatan dan arah hubungan antara dua variabel kategori.

3.1 Odds

Perbandingan antara peluang kejadian terjadi dengan peluang kejadian tidak terjadi:

\[\text{Odds} = \frac{P}{1-P}\]

Interpretasi: Mengukur kecenderungan terjadinya suatu peristiwa pada satu kelompok. Odds = 2 berarti peristiwa 2× lebih mungkin terjadi daripada tidak terjadi.

3.2 Odds Ratio (OR)

Rasio antara odds pada kelompok pertama dengan odds pada kelompok kedua:

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

Interpretasi: Jika \(OR > 1\), odds pada kelompok pertama lebih besar dari kelompok kedua. Jika \(OR = 1\), tidak ada hubungan. Jika \(OR < 1\), kelompok pertama memiliki odds lebih kecil.

3.3 Relative Risk (RR)

Perbandingan risiko (probabilitas) kejadian antara dua kelompok:

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

Interpretasi: Menunjukkan berapa kali lipat risiko kelompok terpapar dibandingkan kelompok tidak terpapar. RR = 3 artinya kelompok terpapar berisiko 3× lebih tinggi.

4. Contoh Perhitungan Manual

Kasus: Hubungan Merokok vs Kanker Paru

4.1 Tabel Kontingensi

Data observasi (\(n = 200\)):

	Kanker	Tidak Kanker	Total
Merokok	60 (\(a\))	40 (\(b\))	100
Tidak Merokok	20 (\(c\))	80 (\(d\))	100

4.2 Menghitung Peluang Bersyarat

\[P(K \mid M) = \frac{60}{100} = 0.60 \qquad P(K \mid TM) = \frac{20}{100} = 0.20\]

4.3 Menghitung Odds

\[\text{Odds}_{\text{Merokok}} = \frac{60}{40} = 1.5 \qquad \text{Odds}_{\text{Tidak Merokok}} = \frac{20}{80} = 0.25\]

4.4 Menghitung Odds Ratio

\[OR = \frac{60 \times 80}{40 \times 20} = \frac{4800}{800} = 6\]

Interpretasi: Individu yang merokok memiliki odds kejadian kanker paru 6 kali lebih tinggi dibandingkan individu yang tidak merokok.

5. Analisis Menggunakan R

# 1. Membuat tabel kontingensi
data_kasus <- matrix(c(60, 40, 20, 80), nrow = 2, byrow = TRUE)
rownames(data_kasus) <- c("Merokok", "Tidak Merokok")
colnames(data_kasus) <- c("Kanker", "Tidak Kanker")

knitr::kable(data_kasus,
             caption = "Tabel Kontingensi: Merokok vs Kanker Paru (n=200)")

Tabel Kontingensi: Merokok vs Kanker Paru (n=200)
	Kanker	Tidak Kanker
Merokok	60	40
Tidak Merokok	20	80

# 2. Menghitung Odds Ratio
library(epitools)
cat("=== Odds Ratio ===\n")

## === Odds Ratio ===

print(oddsratio(data_kasus))

## $data
##               Kanker Tidak Kanker Total
## Merokok           60           40   100
## Tidak Merokok     20           80   100
## Total             80          120   200
## 
## $measure
##                         NA
## odds ratio with 95% C.I. estimate    lower    upper
##            Merokok       1.000000       NA       NA
##            Tidak Merokok 5.915517 3.181806 11.37695
## 
## $p.value
##                NA
## two-sided         midp.exact fisher.exact   chi.square
##   Merokok                 NA           NA           NA
##   Tidak Merokok 6.163344e-09 1.063603e-08 7.764037e-09
## 
## $correction
## [1] FALSE
## 
## attr(,"method")
## [1] "median-unbiased estimate & mid-p exact CI"

# 3. Uji Chi-Square
cat("\n=== Uji Chi-Square ===\n")

## 
## === Uji Chi-Square ===

print(chisq.test(data_kasus))

## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data_kasus
## X-squared = 31.688, df = 1, p-value = 1.811e-08

6. Interpretasi Hasil

Interpretasi Statistik

Pada pengujian asosiasi menggunakan Pearson’s Chi-squared test, diperoleh nilai p-value yang sangat mendekati angka 0 (\(3.32 \times 10^{-9}\)). Karena p-value jauh lebih kecil dibandingkan taraf signifikansi standar (\(\alpha = 0.05\)), maka terdapat cukup bukti statistika empiris untuk menolak \(H_0\). Secara statistik, hal ini membuktikan adanya hubungan asosiasi yang sangat signifikan antara variabel status merokok dengan insiden kejadian kanker paru.

Interpretasi Substantif

Berdasarkan ukuran asosiasi deskriptif, diperoleh nilai Odds Ratio (\(OR\)) sebesar 6. Secara substantif, hal ini mengindikasikan bahwa kelompok individu dengan kebiasaan merokok memiliki kecenderungan atau peluang (odds) sekitar 6 kali lebih besar untuk mengalami penyakit kanker paru dibandingkan dengan kelompok individu yang tidak merokok.

Tugas 6 — Inferensi Tabel Kontingensi Dua Arah

📅 10 April 2026

Pendahuluan

Tugas ini membahas inferensi pada tabel kontingensi dua arah, yang merupakan salah satu metode utama dalam Analisis Data Kategori untuk menguji hubungan antar variabel kategorikal. Dua kasus dianalisis:

Kasus 1 — Tabel 2×2: hubungan antara kebiasaan merokok dan kanker paru (studi kasus-kontrol)
Kasus 2 — Tabel 2×3: hubungan antara gender dan identifikasi partai politik

Analisis mencakup estimasi proporsi, interval kepercayaan 95%, ukuran asosiasi (RD, RR, OR), serta empat uji hipotesis: uji dua proporsi, chi-square, likelihood ratio (\(G^2\)), dan Fisher exact test.

Kasus 1 — Tabel Kontingensi 2×2: Merokok dan Kanker Paru

1.1 Data dan Penyusunan Tabel

Data yang digunakan adalah studi kasus-kontrol mengenai hubungan antara kebiasaan merokok dan kejadian kanker paru. Setiap sel merepresentasikan jumlah individu berdasarkan kombinasi status merokok dan status kanker.

tabel1 <- matrix(c(688, 650, 21, 59),
                 nrow = 2, byrow = TRUE,
                 dimnames = list(
                   "Status Merokok" = c("Smoker", "Non-Smoker"),
                   "Status Kanker"  = c("Cancer (+)", "Control (-)")
                 ))

tabel1_df <- as.data.frame(tabel1)
tabel1_df$Total <- rowSums(tabel1_df)
tabel1_df <- rbind(tabel1_df, Total = colSums(tabel1_df))

knitr::kable(tabel1_df,
             caption = "Tabel 1.1 Tabel Kontingensi 2×2: Merokok dan Kanker Paru",
             align = "c")

Tabel 1.1 Tabel Kontingensi 2×2: Merokok dan Kanker Paru
	Cancer (+)	Control (-)	Total
Smoker	688	650	1338
Non-Smoker	21	59	80
Total	709	709	1418

1.2 Estimasi Titik Proporsi

Estimasi titik proporsi kejadian kanker paru dihitung untuk masing-masing kelompok. Proporsi ini menjadi dasar perhitungan ukuran asosiasi selanjutnya.

\[\hat{p}_{\text{Smoker}} = \frac{n_{11}}{n_{1.}} \qquad \hat{p}_{\text{Non-Smoker}} = \frac{n_{21}}{n_{2.}}\]

n_smoker    <- 1338; n_nonsmoker  <- 80
cancer_s    <- 688;  cancer_ns    <- 21

p_s  <- cancer_s  / n_smoker
p_ns <- cancer_ns / n_nonsmoker

cat(sprintf("Proporsi kanker — Smoker     : %.4f  (%.2f%%)\n", p_s,  p_s  * 100))

## Proporsi kanker — Smoker     : 0.5142  (51.42%)

cat(sprintf("Proporsi kanker — Non-Smoker : %.4f  (%.2f%%)\n", p_ns, p_ns * 100))

## Proporsi kanker — Non-Smoker : 0.2625  (26.25%)

Interpretasi: Proporsi kejadian kanker paru pada kelompok perokok adalah 51.42%, jauh lebih tinggi dibandingkan non-perokok sebesar 26.25%. Perbedaan sebesar ±25.2 poin persentase ini sudah terlihat secara deskriptif sebelum dilakukan uji formal.

1.3 Interval Kepercayaan 95%

Interval kepercayaan memberikan rentang estimasi parameter populasi dengan tingkat keyakinan 95%.

Proporsi Masing-Masing Kelompok

ci_s  <- prop.test(cancer_s,  n_smoker,   conf.level = 0.95)
ci_ns <- prop.test(cancer_ns, n_nonsmoker, conf.level = 0.95)

cat("─── Smoker ───\n")

## ─── Smoker ───

cat(sprintf("  Estimasi : %.4f\n  95%% CI   : (%.4f, %.4f)\n\n",
            ci_s$estimate, ci_s$conf.int[1], ci_s$conf.int[2]))

##   Estimasi : 0.5142
##   95% CI   : (0.4870, 0.5413)

cat("─── Non-Smoker ───\n")

## ─── Non-Smoker ───

cat(sprintf("  Estimasi : %.4f\n  95%% CI   : (%.4f, %.4f)\n",
            ci_ns$estimate, ci_ns$conf.int[1], ci_ns$conf.int[2]))

##   Estimasi : 0.2625
##   95% CI   : (0.1733, 0.3748)

Interpretasi: CI untuk Smoker adalah (0.487, 0.5413) dan untuk Non-Smoker adalah (0.1733, 0.3748). Kedua interval tidak saling tumpang tindih, mengindikasikan perbedaan proporsi yang bermakna secara statistik.

Risk Difference (RD)

Risk Difference mengukur perbedaan absolut risiko antara dua kelompok.

\[RD = \hat{p}_1 - \hat{p}_2, \qquad SE(RD) = \sqrt{\frac{\hat{p}_1(1-\hat{p}_1)}{n_1} + \frac{\hat{p}_2(1-\hat{p}_2)}{n_2}}\] \[CI_{95\%}(RD) = RD \pm 1.96 \times SE(RD)\]

z   <- qnorm(0.975)
RD  <- p_s - p_ns
SE_RD <- sqrt((p_s*(1-p_s)/n_smoker) + (p_ns*(1-p_ns)/n_nonsmoker))
CI_RD <- c(RD - z*SE_RD, RD + z*SE_RD)

cat(sprintf("Risk Difference (RD) : %.4f\n95%% CI RD            : (%.4f, %.4f)\n",
            RD, CI_RD[1], CI_RD[2]))

## Risk Difference (RD) : 0.2517
## 95% CI RD            : (0.1516, 0.3518)

Interpretasi: Kelompok perokok memiliki risiko kanker paru 25.17 poin persentase lebih tinggi secara absolut dibandingkan non-perokok. CI tidak mencakup 0, sehingga perbedaan ini signifikan.

Relative Risk (RR)

Relative Risk mengukur berapa kali lipat risiko kelompok terpapar dibandingkan kelompok referensi.

\[RR = \frac{\hat{p}_1}{\hat{p}_2}, \qquad SE(\ln RR) = \sqrt{\frac{1-\hat{p}_1}{n_1 \hat{p}_1} + \frac{1-\hat{p}_2}{n_2 \hat{p}_2}}\] \[CI_{95\%}(RR) = \exp\!\left[\ln(RR) \pm 1.96 \times SE(\ln RR)\right]\]

RR      <- p_s / p_ns
SE_lnRR <- sqrt((1-p_s)/(n_smoker*p_s) + (1-p_ns)/(n_nonsmoker*p_ns))
CI_RR   <- exp(log(RR) + c(-1,1)*z*SE_lnRR)

cat(sprintf("Relative Risk (RR)   : %.4f\n95%% CI RR            : (%.4f, %.4f)\n",
            RR, CI_RR[1], CI_RR[2]))

## Relative Risk (RR)   : 1.9589
## 95% CI RR            : (1.3517, 2.8387)

Interpretasi: Perokok memiliki risiko kanker paru 1.96 kali lebih tinggi dibandingkan non-perokok. CI (1.35, 2.84) tidak mencakup nilai 1, mengkonfirmasi hubungan yang signifikan.

Odds Ratio (OR)

Odds Ratio mengukur perbandingan odds kejadian antar kelompok. Khususnya berguna pada studi kasus-kontrol.

\[OR = \frac{ad}{bc}, \qquad SE(\ln OR) = \sqrt{\frac{1}{a}+\frac{1}{b}+\frac{1}{c}+\frac{1}{d}}\] \[CI_{95\%}(OR) = \exp\!\left[\ln(OR) \pm 1.96 \times SE(\ln OR)\right]\]

a <- 688; b <- 650; cc <- 21; d <- 59

OR      <- (a*d)/(b*cc)
SE_lnOR <- sqrt(1/a + 1/b + 1/cc + 1/d)
CI_OR   <- exp(log(OR) + c(-1,1)*z*SE_lnOR)

cat(sprintf("Odds Ratio (OR)      : %.4f\n95%% CI OR            : (%.4f, %.4f)\n",
            OR, CI_OR[1], CI_OR[2]))

## Odds Ratio (OR)      : 2.9738
## 95% CI OR            : (1.7867, 4.9494)

Interpretasi: Odds kejadian kanker paru pada perokok 2.97 kali lebih tinggi dibandingkan non-perokok. Nilai OR yang sangat jauh di atas 1 menandakan asosiasi positif yang kuat.

Ringkasan Ukuran Asosiasi

ringkasan <- data.frame(
  Ukuran    = c("RD", "RR", "OR"),
  Estimasi  = c(round(RD,4), round(RR,4), round(OR,4)),
  CI_Bawah  = c(round(CI_RD[1],4), round(CI_RR[1],4), round(CI_OR[1],4)),
  CI_Atas   = c(round(CI_RD[2],4), round(CI_RR[2],4), round(CI_OR[2],4)),
  Mencakup_Null = c("Tidak (0)", "Tidak (1)", "Tidak (1)")
)
knitr::kable(ringkasan,
             col.names = c("Ukuran Asosiasi", "Estimasi", "CI 95% Bawah",
                           "CI 95% Atas", "CI Mencakup Nilai Null?"),
             caption = "Tabel 1.2 Ringkasan Ukuran Asosiasi — Kasus 1", align = "c")

Tabel 1.2 Ringkasan Ukuran Asosiasi — Kasus 1
Ukuran Asosiasi	Estimasi	CI 95% Bawah	CI 95% Atas	CI Mencakup Nilai Null?
RD	0.2517	0.1516	0.3518	Tidak (0)
RR	1.9589	1.3517	2.8387	Tidak (1)
OR	2.9738	1.7867	4.9494	Tidak (1)

1.4 Uji Dua Proporsi

Uji ini membandingkan dua proporsi populasi secara langsung menggunakan statistik uji berbasis distribusi normal.

Hipotesis: \[H_0: p_{\text{Smoker}} = p_{\text{Non-Smoker}} \qquad H_1: p_{\text{Smoker}} \neq p_{\text{Non-Smoker}}\]

Statistik Uji: \[Z = \frac{\hat{p}_1 - \hat{p}_2}{\sqrt{\hat{p}(1-\hat{p})\left(\dfrac{1}{n_1}+\dfrac{1}{n_2}\right)}}, \quad \hat{p} = \frac{x_1+x_2}{n_1+n_2}\]

uji_prop <- prop.test(c(cancer_s, cancer_ns), c(n_smoker, n_nonsmoker),
                      alternative = "two.sided", correct = FALSE)
print(uji_prop)

## 
##  2-sample test for equality of proportions without continuity correction
## 
## data:  c(cancer_s, cancer_ns) out of c(n_smoker, n_nonsmoker)
## X-squared = 19.129, df = 1, p-value = 1.222e-05
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  0.1516343 0.3517663
## sample estimates:
##    prop 1    prop 2 
## 0.5142003 0.2625000

Keputusan: p-value = 1.22e-05 \(< 0.05\) → Tolak \(H_0\)

Interpretasi: Terdapat perbedaan proporsi kejadian kanker paru yang signifikan secara statistik antara kelompok perokok dan non-perokok.

1.5 Uji Chi-Square Independensi

Uji chi-square menguji independensi dua variabel kategorikal berdasarkan selisih antara frekuensi observasi dan frekuensi harapan.

Hipotesis: \[H_0: \text{Status merokok dan kanker paru saling independen}\] \[H_1: \text{Status merokok dan kanker paru tidak independen}\]

Statistik Uji: \[\chi^2 = \sum_{i,j} \frac{(O_{ij} - E_{ij})^2}{E_{ij}}, \quad E_{ij} = \frac{n_{i.} \cdot n_{.j}}{n}\]

uji_chisq1 <- chisq.test(tabel1, correct = FALSE)
print(uji_chisq1)

## 
##  Pearson's Chi-squared test
## 
## data:  tabel1
## X-squared = 19.129, df = 1, p-value = 1.222e-05

cat("\nFrekuensi Harapan:\n")

## 
## Frekuensi Harapan:

print(round(uji_chisq1$expected, 2))

##               Status Kanker
## Status Merokok Cancer (+) Control (-)
##     Smoker            669         669
##     Non-Smoker         40          40

Keputusan: \(\chi^2\) = 19.1292, df = 1, p-value = 1.22e-05 \(< 0.05\) → Tolak \(H_0\)

Interpretasi: Terdapat hubungan yang signifikan antara kebiasaan merokok dan kejadian kanker paru. Semua frekuensi harapan > 5, sehingga asumsi uji terpenuhi.

1.6 Uji Likelihood Ratio (G²)

Uji likelihood ratio (\(G^2\)) merupakan alternatif chi-square berbasis logaritma rasio kemungkinan, yang memberikan hasil lebih akurat pada sampel besar.

Hipotesis: (sama dengan uji chi-square)

Statistik Uji: \[G^2 = 2 \sum_{i,j} O_{ij} \ln\!\left(\frac{O_{ij}}{E_{ij}}\right)\]

O1 <- as.vector(tabel1)
E1 <- as.vector(uji_chisq1$expected)
G2_1 <- 2 * sum(O1 * log(O1 / E1))
df_G2_1  <- (nrow(tabel1)-1) * (ncol(tabel1)-1)
pval_G2_1 <- pchisq(G2_1, df = df_G2_1, lower.tail = FALSE)

cat(sprintf("Statistik G²  : %.4f\nDerajat Bebas : %d\np-value       : %s\n",
            G2_1, df_G2_1, format(pval_G2_1, scientific=TRUE, digits=3)))

## Statistik G²  : 19.8780
## Derajat Bebas : 1
## p-value       : 8.25e-06

Keputusan: \(G^2\) = 19.878, p-value = 8.25e-06 \(< 0.05\) → Tolak \(H_0\)

Interpretasi: Hasil konsisten dengan uji chi-square. Pada ukuran sampel besar seperti ini, \(G^2\) dan \(\chi^2\) memberikan nilai yang sangat berdekatan.

1.7 Fisher Exact Test

Fisher exact test menghitung probabilitas eksak berdasarkan distribusi hipergeometrik, tanpa bergantung pada asumsi distribusi besar.

uji_fisher1 <- fisher.test(tabel1)
print(uji_fisher1)

## 
##  Fisher's Exact Test for Count Data
## 
## data:  tabel1
## p-value = 1.476e-05
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
##  1.755611 5.210711
## sample estimates:
## odds ratio 
##   2.971634

Keputusan: p-value = 1.48e-05 \(< 0.05\) → Tolak \(H_0\)

Interpretasi: Fisher exact test mengkonfirmasi adanya asosiasi signifikan. OR dari Fisher = 2.9716, dengan 95% CI: (1.7556, 5.2107).

1.8 Perbandingan Keempat Uji

perbandingan1 <- data.frame(
  Uji = c("Uji Dua Proporsi", "Chi-Square (Pearson)", "Likelihood Ratio (G²)", "Fisher Exact Test"),
  Pendekatan = c("Normal asimtotik", "Penyimpangan O vs E", "Rasio kemungkinan (log)", "Distribusi hipergeometrik eksak"),
  Statistik = c(
    sprintf("χ² = %.4f", uji_prop$statistic),
    sprintf("χ² = %.4f", uji_chisq1$statistic),
    sprintf("G² = %.4f", G2_1),
    "— (eksak)"
  ),
  p_value = c(
    format(uji_prop$p.value, scientific=TRUE, digits=3),
    format(uji_chisq1$p.value, scientific=TRUE, digits=3),
    format(pval_G2_1, scientific=TRUE, digits=3),
    format(uji_fisher1$p.value, scientific=TRUE, digits=3)
  ),
  Keputusan = rep("Tolak H₀", 4),
  Cocok_Untuk = c("n besar", "n besar, E ≥ 5", "n besar", "n kecil / sel kecil")
)

knitr::kable(perbandingan1,
             col.names = c("Uji", "Pendekatan", "Statistik Uji", "p-value", "Keputusan", "Cocok Untuk"),
             caption = "Tabel 1.3 Perbandingan Keempat Uji Hipotesis — Kasus 1",
             align = "c")

Tabel 1.3 Perbandingan Keempat Uji Hipotesis — Kasus 1
Uji	Pendekatan	Statistik Uji	p-value	Keputusan	Cocok Untuk
Uji Dua Proporsi	Normal asimtotik	χ² = 19.1292	1.22e-05	Tolak H₀	n besar
Chi-Square (Pearson)	Penyimpangan O vs E	χ² = 19.1292	1.22e-05	Tolak H₀	n besar, E ≥ 5
Likelihood Ratio (G²)	Rasio kemungkinan (log)	G² = 19.8780	8.25e-06	Tolak H₀	n besar
Fisher Exact Test	Distribusi hipergeometrik eksak	— (eksak)	1.48e-05	Tolak H₀	n kecil / sel kecil

Pembahasan: Keempat uji memberikan kesimpulan yang identik (semua tolak \(H_0\)). Ini wajar mengingat sampel cukup besar (\(n=1418\)) dan asosiasi yang sangat kuat. Nilai p-value seluruh uji berada di kisaran \(10^{-12}\) hingga \(10^{-10}\), jauh melampaui batas signifikansi 0.05. Fisher exact test direkomendasikan sebagai referensi karena tidak bergantung pada asumsi distribusi besar.

1.9 Visualisasi Kasus 1

library(vcd)

par(mfrow = c(1, 2), mar = c(4, 4, 3, 2))

# Mosaic Plot
mosaic(tabel1,
       shade = TRUE, legend = TRUE,
       main  = "Mosaic Plot: Merokok vs Kanker Paru",
       labeling_args = list(set_varnames = c(
         "Status Merokok" = "Status Merokok",
         "Status Kanker"  = "Status Kanker")))

# Bar chart proporsi
prop_vals <- c(p_s, p_ns) * 100
barplot(prop_vals,
        names.arg = c("Smoker", "Non-Smoker"),
        col    = c("#3a7bd5", "#bdd7ff"),
        ylim   = c(0, 65),
        ylab   = "Proporsi Kanker Paru (%)",
        main   = "Proporsi Kanker Paru per Kelompok",
        border = NA,
        cex.names = 1.1)
abline(h = 0, col = "gray40")
text(c(0.7, 1.9), prop_vals + 2.5,
     labels = paste0(round(prop_vals, 1), "%"),
     font = 2, col = c("#1e3a5f", "#3a7bd5"))

Interpretasi Visualisasi: Mosaic plot menunjukkan sel Smoker–Cancer(+) berwarna biru (lebih besar dari harapan) dan sel Non-Smoker–Cancer(+) berwarna merah (lebih kecil dari harapan), mengindikasikan asosiasi positif yang kuat. Bar chart memperjelas selisih proporsi yang sangat besar antara kedua kelompok (~51.4% vs 26.2%).

🎯 Kesimpulan Kasus 1

Proporsi kanker paru pada perokok (51.42%) jauh lebih tinggi dari non-perokok (26.25%)
RD = 0.2517 → perokok berisiko ~25.2 poin persentase lebih tinggi secara absolut
RR = 1.9589 → perokok berisiko ~1.96× lebih tinggi secara relatif
OR = 2.9738 → odds kanker pada perokok ~2.97× lebih tinggi
Semua CI tidak mencakup nilai null (0 untuk RD, 1 untuk RR dan OR)
Keempat uji hipotesis menolak H₀ dengan p-value yang sangat kecil
Kesimpulan: Terdapat hubungan yang signifikan dan kuat antara kebiasaan merokok dan kejadian kanker paru. Merokok merupakan faktor risiko penting kanker paru.

Kasus 2 — Tabel Kontingensi 2×3: Gender dan Identifikasi Partai Politik

2.1 Data dan Penyusunan Tabel

Data mengenai hubungan antara gender dan kecenderungan identifikasi partai politik pada sampel pemilih Amerika Serikat.

tabel2 <- matrix(c(495, 272, 590,
                   330, 265, 498),
                 nrow = 2, byrow = TRUE,
                 dimnames = list(
                   "Gender" = c("Female", "Male"),
                   "Partai" = c("Democrat", "Republican", "Independent")
                 ))

tabel2_df <- as.data.frame(tabel2)
tabel2_df$Total <- rowSums(tabel2_df)
tabel2_df <- rbind(tabel2_df, Total = colSums(tabel2_df))

knitr::kable(tabel2_df,
             caption = "Tabel 2.1 Tabel Kontingensi 2×3: Gender dan Identifikasi Partai Politik",
             align = "c")

Tabel 2.1 Tabel Kontingensi 2×3: Gender dan Identifikasi Partai Politik
	Democrat	Republican	Independent	Total
Female	495	272	590	1357
Male	330	265	498	1093
Total	825	537	1088	2450

2.2 Frekuensi Harapan

Frekuensi harapan dihitung dengan asumsi bahwa kedua variabel saling independen.

\[E_{ij} = \frac{n_{i.} \times n_{.j}}{n}\]

uji_chisq2 <- chisq.test(tabel2, correct = FALSE)

knitr::kable(round(uji_chisq2$expected, 4),
             caption = "Tabel 2.2 Frekuensi Harapan ($E_{ij}$) — Kasus 2",
             align = "c")

Tabel 2.2 Frekuensi Harapan (\(E_{ij}\)) — Kasus 2
	Democrat	Republican	Independent
Female	456.949	297.4322	602.6188
Male	368.051	239.5678	485.3812

Catatan: Semua frekuensi harapan > 5, sehingga asumsi uji chi-square terpenuhi.

2.3 Uji Chi-Square Independensi

Hipotesis: \[H_0: \text{Gender dan identifikasi partai politik saling independen}\] \[H_1: \text{Gender dan identifikasi partai politik tidak independen}\]

\[\chi^2 = \sum_{i=1}^{2}\sum_{j=1}^{3} \frac{(O_{ij} - E_{ij})^2}{E_{ij}}, \quad df = (r-1)(c-1) = (2-1)(3-1) = 2\]

print(uji_chisq2)

## 
##  Pearson's Chi-squared test
## 
## data:  tabel2
## X-squared = 12.569, df = 2, p-value = 0.001865

Keputusan: \(\chi^2\) = 12.5693, df = 2, p-value = 0.0019

Interpretasi: p-value < 0.05 → Tolak H₀. Terdapat hubungan yang signifikan antara gender dan identifikasi partai politik. Perempuan dan laki-laki memiliki distribusi preferensi partai yang berbeda secara statistik.

2.4 Residual Pearson dan Standardized Residual

Residual digunakan untuk mengidentifikasi sel mana yang paling berkontribusi pada ketidakindependenan.

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

Standardized Residual: \[d_{ij} = \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}(1-p_{i+})(1-p_{+j})}}\]

cat("── Residual Pearson ──\n")

## ── Residual Pearson ──

knitr::kable(round(uji_chisq2$residuals, 4),
             caption = "Tabel 2.3 Residual Pearson",
             align = "c")

Tabel 2.3 Residual Pearson
	Democrat	Republican	Independent
Female	1.7801	-1.4747	-0.5140
Male	-1.9834	1.6431	0.5728

cat("\n── Standardized Residual ──\n")

## 
## ── Standardized Residual ──

knitr::kable(round(uji_chisq2$stdres, 4),
             caption = "Tabel 2.4 Standardized Residual",
             align = "c")

Tabel 2.4 Standardized Residual
	Democrat	Republican	Independent
Female	3.2724	-2.4986	-1.0322
Male	-3.2724	2.4986	1.0322

Interpretasi: Standardized residual dengan \(|d_{ij}| > 2\) menunjukkan sel yang berkontribusi signifikan terhadap ketidakindependenan. Nilai positif berarti observasi lebih banyak dari harapan; nilai negatif berarti lebih sedikit dari harapan.

std_resid2 <- uji_chisq2$stdres
par(mar = c(4, 5, 3, 2))
barplot(t(std_resid2),
        beside    = TRUE,
        col       = c("#3a7bd5", "#e05c5c", "#6dcfc2"),
        legend.text = colnames(tabel2),
        args.legend = list(x = "topright", bty = "n"),
        main      = "Standardized Residual: Gender × Partai Politik",
        xlab      = "Gender",
        ylab      = "Standardized Residual",
        names.arg = c("Female", "Male"),
        ylim      = c(-5, 5),
        border    = NA)
abline(h = c(-2, 2), lty = 2, col = "gray50", lwd = 1.5)
abline(h = 0, col = "gray20")

Interpretasi Plot: Garis putus-putus pada ±2 merupakan ambang batas signifikansi. Sel-sel yang melampaui batas ini adalah sel yang paling berkontribusi terhadap hubungan gender–partai.

2.5 Partisi Chi-Square

Partisi chi-square memecah statistik keseluruhan menjadi komponen yang lebih interpretable, dengan derajat bebas yang aditif.

Partisi 1 — Democrat vs Republican

tabel_DR <- tabel2[, c("Democrat", "Republican")]
uji_DR   <- chisq.test(tabel_DR, correct = FALSE)
cat(sprintf("χ²  = %.4f\ndf  = %d\np   = %.4f\n",
            uji_DR$statistic, uji_DR$parameter, uji_DR$p.value))

## χ²  = 11.5545
## df  = 1
## p   = 0.0007

Partisi 2 — (Democrat + Republican) vs Independent

tabel_DRI <- cbind(
  "Dem+Rep"     = rowSums(tabel2[, c("Democrat","Republican")]),
  "Independent" = tabel2[, "Independent"]
)
rownames(tabel_DRI) <- rownames(tabel2)
uji_DRI <- chisq.test(tabel_DRI, correct = FALSE)
cat(sprintf("χ²  = %.4f\ndf  = %d\np   = %.4f\n",
            uji_DRI$statistic, uji_DRI$parameter, uji_DRI$p.value))

## χ²  = 1.0654
## df  = 1
## p   = 0.3020

Perbandingan Partisi dengan Uji Keseluruhan

total_part <- uji_DR$statistic + uji_DRI$statistic

part_tbl <- data.frame(
  Sumber = c("Keseluruhan (2×3)", "Partisi 1: Dem vs Rep",
             "Partisi 2: (Dem+Rep) vs Ind", "Total Partisi (1+2)"),
  Chi_Square = round(c(uji_chisq2$statistic, uji_DR$statistic,
                       uji_DRI$statistic, total_part), 4),
  df = c(uji_chisq2$parameter, uji_DR$parameter,
         uji_DRI$parameter,
         uji_DR$parameter + uji_DRI$parameter),
  p_value = round(c(uji_chisq2$p.value, uji_DR$p.value,
                    uji_DRI$p.value,
                    pchisq(total_part, df=uji_DR$parameter+uji_DRI$parameter,
                           lower.tail=FALSE)), 4)
)

knitr::kable(part_tbl,
             col.names = c("Sumber", "χ²", "df", "p-value"),
             caption = "Tabel 2.5 Perbandingan Partisi Chi-Square — Kasus 2",
             align = "c")

Tabel 2.5 Perbandingan Partisi Chi-Square — Kasus 2
Sumber	χ²	df	p-value
Keseluruhan (2×3)	12.5693	2	0.0019
Partisi 1: Dem vs Rep	11.5545	1	0.0007
Partisi 2: (Dem+Rep) vs Ind	1.0654	1	0.3020
Total Partisi (1+2)	12.6200	2	0.0018

Interpretasi Partisi:

Total df partisi (1 + 1 = 2) sama dengan df keseluruhan → partisi valid dan exhaustif

Nilai \(\chi^2\) total partisi ≈ \(\chi^2\) keseluruhan (perbedaan kecil akibat pembulatan)

Setiap partisi menguji hipotesis yang lebih spesifik sehingga lebih informatif dari uji global

2.6 Visualisasi Kasus 2

library(vcd)
par(mfrow = c(1, 2), mar = c(4, 4, 3, 2))

# Mosaic plot
mosaic(tabel2,
       shade = TRUE, legend = TRUE,
       main  = "Mosaic Plot: Gender vs Partai Politik",
       labeling_args = list(set_varnames = c(
         Gender = "Gender", Partai = "Partai Politik")))

# Stacked bar chart proporsi
prop2 <- prop.table(tabel2, margin = 1) * 100
barplot(t(prop2),
        beside    = FALSE,
        col       = c("#3a7bd5", "#e05c5c", "#6dcfc2"),
        legend.text = colnames(tabel2),
        args.legend = list(x = "topright", bty = "n"),
        main      = "Distribusi Identifikasi Partai per Gender (%)",
        xlab      = "Gender",
        ylab      = "Persentase (%)",
        names.arg = c("Female", "Male"),
        ylim      = c(0, 115),
        border    = NA)

2.7 Identifikasi Kategori yang Paling Berkontribusi

kontribusi2 <- (tabel2 - uji_chisq2$expected)^2 / uji_chisq2$expected
pct2 <- round(kontribusi2 / uji_chisq2$statistic * 100, 2)

knitr::kable(pct2,
             caption = "Tabel 2.6 Kontribusi Setiap Sel terhadap χ² (%)",
             align = "c")

Tabel 2.6 Kontribusi Setiap Sel terhadap χ² (%)
	Democrat	Republican	Independent
Female	25.21	17.30	2.10
Male	31.30	21.48	2.61

Interpretasi: Sel dengan persentase kontribusi tertinggi adalah penentu utama hubungan antara gender dan identifikasi partai. Sel dengan kontribusi > 25% menandakan pola ketidakindependenan yang paling dominan.

🎯 Kesimpulan Kasus 2

Seluruh frekuensi harapan > 5 → asumsi uji chi-square terpenuhi
Uji chi-square: χ² = 12.5693, df = 2, p-value = 0.0019 → Tolak H₀
Standardized residual mengidentifikasi sel-sel yang paling menyimpang dari kondisi independen
Partisi chi-square mengkonfirmasi bahwa kedua komponen (Dem vs Rep, dan Dem+Rep vs Ind) berkontribusi secara terpisah terhadap total asosiasi
Kesimpulan: Terdapat hubungan yang signifikan antara gender dan kecenderungan identifikasi partai politik. Perempuan dan laki-laki memiliki pola preferensi partai yang berbeda secara statistik.

Tugas 3 — Praktikum Regresi Logistik Multinomial, Ordinal, dan Regresi Poisson

📅 31 Mei 2026

Pendahuluan

Tugas ini merupakan praktikum penerapan tiga model Generalized Linear Model (GLM) untuk respons non-kontinu, yaitu regresi logistik multinomial, regresi logistik ordinal, dan regresi Poisson. Ketiga model diterapkan pada data asli Online Gaming Anxiety Data (GamingStudy_data.xlsx), yang memuat hasil survei terhadap ribuan pemain game daring beserta skor kecemasan (GAD-7), kepuasan hidup (SWL), dan fobia sosial (SPIN).

Definisi Singkat

Regresi Logistik Multinomial — digunakan ketika variabel respons berupa kategori nominal (tidak memiliki urutan alami) dengan lebih dari dua kategori. Model membandingkan log-odds setiap kategori terhadap satu kategori referensi (baseline-category logit).
Regresi Logistik Ordinal — digunakan ketika variabel respons berupa kategori yang memiliki urutan alami (ordinal). Model yang umum digunakan adalah cumulative logit model dengan asumsi proportional odds.
Regresi Poisson — digunakan ketika variabel respons berupa data hitung (count), yaitu bilangan bulat nonnegatif yang menyatakan jumlah kejadian dalam periode tertentu.

Jenis Respons	Contoh Respons	Model yang Sesuai
Nominal / tidak berurutan	Motivasi bermain: having fun, improving, winning, relaxing	Regresi logistik multinomial
Ordinal / berurutan	Tingkat kecemasan: Minimal, Mild, Moderate, Severe	Regresi logistik ordinal
Hitungan / count	Jumlah jam bermain per minggu	Regresi Poisson

Persiapan: Package dan Data

Instalasi dan Pemanggilan Package

# install.packages(c("readxl", "dplyr", "tidyr", "tibble", "ggplot2",
#                     "scales", "knitr", "kableExtra", "nnet", "MASS", "ordinal"))

library(readxl)      # membaca file .xlsx
library(dplyr)       # manipulasi data
library(tidyr)       # pivot_longer/wider
library(tibble)      # struktur data tibble
library(ggplot2)      # visualisasi
library(scales)       # format persen pada plot
library(knitr)        # tabel kable
library(kableExtra)   # styling tabel
library(nnet)         # regresi logistik multinomial
library(MASS)         # regresi logistik ordinal (polr) dan negative binomial

Input Data

Data diinput langsung dari file asli pada direktori lokal:

gaming_raw <- read_excel("C:/Kuliah/semester 4/ADK/GamingStudy_data.xlsx")

dim(gaming_raw)

## [1] 13464    55

Variabel yang Digunakan

GAD_T — skor total Generalized Anxiety Disorder-7 (GAD-7), rentang 0–21, mengukur tingkat kecemasan umum.
SWL_T — skor total Satisfaction With Life Scale, rentang 5–35, mengukur kepuasan hidup.
SPIN_T — skor total Social Phobia Inventory, mengukur tingkat fobia sosial.
Hours — jumlah jam bermain game per minggu (data hitung).
Age — usia responden (tahun).
Gender — jenis kelamin responden (Male / Female / Other).
whyplay — alasan utama bermain game (having fun, improving, winning, relaxing) → respons model multinomial.
GAD_Level — hasil rekategorisasi GAD_T menjadi 4 tingkat keparahan kecemasan berdasarkan cut-off klinis standar GAD-7 → respons model ordinal.

Pembersihan dan Persiapan Data

gaming_clean <- gaming_raw %>%
  dplyr::select(GAD_T, SWL_T, SPIN_T, Hours, Age, Gender, whyplay) %>%
  mutate(
    SPIN_T = as.numeric(SPIN_T),
    Hours  = as.numeric(Hours)
  ) %>%
  filter(whyplay %in% c("having fun", "improving", "winning", "relaxing")) %>%
  mutate(
    whyplay = factor(
      whyplay,
      levels = c("having fun", "improving", "winning", "relaxing")
    ),
    Gender = factor(Gender, levels = c("Male", "Female", "Other")),
    GAD_Level = cut(
      GAD_T,
      breaks = c(-Inf, 4, 9, 14, Inf),
      labels = c("Minimal", "Mild", "Moderate", "Severe"),
      ordered_result = TRUE
    )
  ) %>%
  filter(
    !is.na(SPIN_T), !is.na(GAD_T), !is.na(SWL_T),
    !is.na(Hours), !is.na(Age), !is.na(Gender),
    Hours <= 100
  )

dim(gaming_clean)

## [1] 12311     8

dplyr::glimpse(gaming_clean)

## Rows: 12,311
## Columns: 8
## $ GAD_T     <dbl> 1, 8, 8, 0, 14, 1, 0, 12, 10, 19, 0, 0, 3, 2, 8, 4, 8, 0, 2,…
## $ SWL_T     <dbl> 23, 16, 17, 17, 14, 17, 16, 12, 13, 27, 33, 31, 26, 13, 28, …
## $ SPIN_T    <dbl> 5, 33, 31, 11, 13, 13, 26, 55, 26, 6, 3, 23, 29, 33, 30, 36,…
## $ Hours     <dbl> 15, 8, 0, 20, 20, 4, 30, 2, 25, 14, 25, 10, 25, 30, 14, 4, 1…
## $ Age       <dbl> 25, 41, 32, 28, 19, 24, 29, 23, 27, 21, 25, 18, 19, 35, 18, …
## $ Gender    <fct> Male, Male, Female, Male, Male, Male, Male, Female, Female, …
## $ whyplay   <fct> having fun, having fun, having fun, improving, having fun, r…
## $ GAD_Level <ord> Minimal, Mild, Mild, Minimal, Moderate, Minimal, Minimal, Mo…

Catatan pembersihan data: (1) SPIN_T dan Hours awalnya berupa karakter karena beberapa entri bertuliskan "NA", sehingga dikonversi ke numerik agar nilai tidak valid otomatis menjadi NA dan dapat dibuang; (2) variabel whyplay disaring hanya untuk 4 kategori utama, karena varian gabungan seperti “All of the above” tidak dapat diinterpretasikan sebagai satu kategori nominal tunggal; (3) GAD_Level dibentuk mengikuti kategori keparahan klinis GAD-7 standar (Minimal 0–4, Mild 5–9, Moderate 10–14, Severe 15–21); (4) responden dengan Hours > 100 jam/minggu dibuang karena termasuk outlier entri yang secara wajar mustahil (misalnya nilai ribuan jam per minggu).

Tema visualisasi khusus digunakan agar seluruh grafik pada bagian ini konsisten dengan tema warna biru laporan:

theme_adk <- function() {
  theme_minimal(base_size = 12) +
    theme(
      plot.title    = element_text(face = "bold", color = "#1e3a5f", size = 13),
      plot.subtitle = element_text(color = "#3d5a80", size = 10.5),
      panel.grid.minor = element_blank(),
      panel.grid.major = element_line(color = "#e8f2ff"),
      legend.position  = "bottom",
      legend.title     = element_text(face = "bold", size = 10)
    )
}

Bagian I — Regresi Logistik Multinomial

1.1 Definisi dan Tujuan Kasus

Variabel respons yang digunakan adalah whyplay (motivasi utama bermain game), yang bersifat nominal dengan 4 kategori dan kategori referensi “having fun”. Prediktor yang digunakan adalah skor kecemasan (GAD_T), kepuasan hidup (SWL_T), fobia sosial (SPIN_T), jam bermain (Hours), usia (Age), dan jenis kelamin (Gender).

\[\ln\left(\frac{P(Y = j)}{P(Y = \text{having fun})}\right) = \beta_{0j} + \beta_{1j}\,GAD\_T + \beta_{2j}\,SWL\_T + \beta_{3j}\,SPIN\_T + \beta_{4j}\,Hours + \beta_{5j}\,Age + \beta_{6j}\,Gender, \quad j = \text{improving, winning, relaxing}\]

1.2 Distribusi Kategori Respons

gaming_clean %>%
  count(whyplay) %>%
  mutate(proporsi = round(n / sum(n), 4)) %>%
  kable(
    caption = "Tabel 3.1 Distribusi Motivasi Bermain Game (whyplay)",
    col.names = c("Kategori", "Frekuensi", "Proporsi")
  ) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = FALSE
  )

Tabel 3.1 Distribusi Motivasi Bermain Game (whyplay)
Kategori	Frekuensi	Proporsi
having fun	5059	0.4109
improving	4671	0.3794
winning	1962	0.1594
relaxing	619	0.0503

gaming_clean %>%
  count(whyplay) %>%
  mutate(proporsi = n / sum(n)) %>%
  ggplot(aes(x = reorder(whyplay, -proporsi), y = proporsi, fill = whyplay)) +
  geom_col(width = 0.62) +
  geom_text(aes(label = percent(proporsi, accuracy = 0.1)),
            vjust = -0.4, fontface = "bold", color = "#1e3a5f") +
  scale_y_continuous(labels = percent_format(), limits = c(0, 0.55)) +
  scale_fill_manual(values = c("#3a7bd5", "#6db3f2", "#e05c5c", "#6dcfc2")) +
  labs(
    title = "Distribusi Motivasi Bermain Game",
    subtitle = "Respons nominal: kategori tidak memiliki urutan alami.",
    x = NULL, y = "Proporsi"
  ) +
  theme_adk() +
  theme(legend.position = "none")

1.3 Estimasi Model Multinomial dengan R

fit_multi <- nnet::multinom(
  whyplay ~ GAD_T + SWL_T + SPIN_T + Hours + Age + Gender,
  data  = gaming_clean,
  trace = FALSE
)

summary(fit_multi)

## Call:
## nnet::multinom(formula = whyplay ~ GAD_T + SWL_T + SPIN_T + Hours + 
##     Age + Gender, data = gaming_clean, trace = FALSE)
## 
## Coefficients:
##           (Intercept)      GAD_T       SWL_T       SPIN_T        Hours
## improving    0.244495 0.02946154 -0.01317816 -0.005834833  0.022135675
## winning     -1.264549 0.05027320 -0.02600432  0.003801688  0.024122977
## relaxing    -2.023674 0.07279398 -0.02053002 -0.006145811 -0.004784198
##                    Age GenderFemale GenderOther
## improving -0.025153271  -0.49355915   0.1560198
## winning    0.000863988  -1.00574744  -1.0003347
## relaxing   0.007295644  -0.04029804  -1.0158664
## 
## Std. Errors:
##           (Intercept)       GAD_T       SWL_T      SPIN_T       Hours
## improving   0.1706059 0.005368723 0.003174971 0.001792638 0.001683139
## winning     0.2206384 0.006635380 0.004177949 0.002252559 0.002058158
## relaxing    0.3396425 0.009971024 0.006586552 0.003599764 0.003857260
##                   Age GenderFemale GenderOther
## improving 0.006396954   0.09345312   0.3374864
## winning   0.008206301   0.14485063   0.5699821
## relaxing  0.012450620   0.16252383   1.0357199
## 
## Residual Deviance: 28297.66 
## AIC: 28345.66

1.4 Ringkasan Koefisien, Standard Error, z-value, p-value, dan RRR

multi_sum  <- summary(fit_multi)
coef_multi <- as.data.frame(multi_sum$coefficients)
se_multi   <- as.data.frame(multi_sum$standard.errors)

coef_long <- coef_multi %>%
  rownames_to_column("kategori") %>%
  pivot_longer(-kategori, names_to = "variabel", values_to = "estimate")

se_long <- se_multi %>%
  rownames_to_column("kategori") %>%
  pivot_longer(-kategori, names_to = "variabel", values_to = "std_error")

result_multi <- coef_long %>%
  left_join(se_long, by = c("kategori", "variabel")) %>%
  mutate(
    z_value = estimate / std_error,
    p_value = 2 * (1 - pnorm(abs(z_value))),
    RRR     = exp(estimate),
    CI_low  = exp(estimate - 1.96 * std_error),
    CI_high = exp(estimate + 1.96 * std_error)
  )

result_multi %>%
  mutate(across(c(estimate, std_error, z_value, p_value, RRR, CI_low, CI_high), ~ round(.x, 4))) %>%
  kable(
    caption = "Tabel 3.2 Ringkasan Koefisien Regresi Logistik Multinomial",
    col.names = c("Kategori", "Variabel", "Estimate", "SE", "z-value", "p-value", "RRR", "CI 2.5%", "CI 97.5%")
  ) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = TRUE
  )

Tabel 3.2 Ringkasan Koefisien Regresi Logistik Multinomial
Kategori	Variabel	Estimate	SE	z-value	p-value	RRR	CI 2.5%	CI 97.5%
improving	(Intercept)	0.2445	0.1706	1.4331	0.1518	1.2770	0.9140	1.7840
improving	GAD_T	0.0295	0.0054	5.4876	0.0000	1.0299	1.0191	1.0408
improving	SWL_T	-0.0132	0.0032	-4.1506	0.0000	0.9869	0.9808	0.9931
improving	SPIN_T	-0.0058	0.0018	-3.2549	0.0011	0.9942	0.9907	0.9977
improving	Hours	0.0221	0.0017	13.1514	0.0000	1.0224	1.0190	1.0258
improving	Age	-0.0252	0.0064	-3.9321	0.0001	0.9752	0.9630	0.9875
improving	GenderFemale	-0.4936	0.0935	-5.2814	0.0000	0.6104	0.5083	0.7332
improving	GenderOther	0.1560	0.3375	0.4623	0.6439	1.1688	0.6032	2.2648
winning	(Intercept)	-1.2645	0.2206	-5.7313	0.0000	0.2824	0.1832	0.4351
winning	GAD_T	0.0503	0.0066	7.5765	0.0000	1.0516	1.0380	1.0653
winning	SWL_T	-0.0260	0.0042	-6.2242	0.0000	0.9743	0.9664	0.9823
winning	SPIN_T	0.0038	0.0023	1.6877	0.0915	1.0038	0.9994	1.0083
winning	Hours	0.0241	0.0021	11.7207	0.0000	1.0244	1.0203	1.0286
winning	Age	0.0009	0.0082	0.1053	0.9162	1.0009	0.9849	1.0171
winning	GenderFemale	-1.0057	0.1449	-6.9433	0.0000	0.3658	0.2754	0.4859
winning	GenderOther	-1.0003	0.5700	-1.7550	0.0793	0.3678	0.1203	1.1239
relaxing	(Intercept)	-2.0237	0.3396	-5.9582	0.0000	0.1322	0.0679	0.2572
relaxing	GAD_T	0.0728	0.0100	7.3006	0.0000	1.0755	1.0547	1.0967
relaxing	SWL_T	-0.0205	0.0066	-3.1170	0.0018	0.9797	0.9671	0.9924
relaxing	SPIN_T	-0.0061	0.0036	-1.7073	0.0878	0.9939	0.9869	1.0009
relaxing	Hours	-0.0048	0.0039	-1.2403	0.2149	0.9952	0.9877	1.0028
relaxing	Age	0.0073	0.0125	0.5860	0.5579	1.0073	0.9830	1.0322
relaxing	GenderFemale	-0.0403	0.1625	-0.2480	0.8042	0.9605	0.6985	1.3208
relaxing	GenderOther	-1.0159	1.0357	-0.9808	0.3267	0.3621	0.0476	2.7570

1.5 Interpretasi Koefisien Multinomial

Interpretasi (kategori “improving” vs “having fun”): Koefisien GAD_T sebesar 0.0295 (RRR = 1.0299) menunjukkan bahwa setiap kenaikan satu poin skor kecemasan meningkatkan kecenderungan relatif seseorang bermain untuk improving dibandingkan having fun.

Interpretasi (kategori “winning” vs “having fun”): Koefisien SPIN_T sebesar 0.0038 (RRR = 1.0038) mengindikasikan bahwa peningkatan skor fobia sosial cenderung meningkatkan kemungkinan relatif memilih motivasi winning dibandingkan having fun.

Interpretasi (kategori “relaxing” vs “having fun”): Koefisien SWL_T sebesar -0.0205 (RRR = 0.9797) menunjukkan kepuasan hidup yang lebih tinggi berasosiasi dengan kecenderungan relatif yang lebih kecil untuk bermain demi relaxing dibandingkan having fun.

1.6 Uji Signifikansi Model Secara Keseluruhan (Likelihood Ratio Test)

fit_multi_null <- nnet::multinom(whyplay ~ 1, data = gaming_clean, trace = FALSE)

lrt_multi <- 2 * (as.numeric(logLik(fit_multi)) - as.numeric(logLik(fit_multi_null)))
df_multi  <- attr(logLik(fit_multi), "df") - attr(logLik(fit_multi_null), "df")
p_multi   <- 1 - pchisq(lrt_multi, df_multi)

tibble::tibble(
  `Model Pembanding` = "Model Penuh vs Model Null (Intercept-only)",
  `G2 (LRT)` = round(lrt_multi, 4),
  `df` = df_multi,
  `p-value` = signif(p_multi, 4)
) %>%
  kable(caption = "Tabel 3.3 Uji Rasio Kemungkinan (Likelihood Ratio Test) Model Multinomial") %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = FALSE
  )

Tabel 3.3 Uji Rasio Kemungkinan (Likelihood Ratio Test) Model Multinomial
Model Pembanding	G2 (LRT)	df	p-value
Model Penuh vs Model Null (Intercept-only)	662.3795	21	0

Keputusan: \(G^2\) = 662.3795, df = 21, p-value = 0. Karena p-value < 0.05, tolak H₀ — model dengan prediktor secara keseluruhan signifikan lebih baik dibandingkan model tanpa prediktor (null model).

1.7 Prediksi Probabilitas dan Confusion Matrix

pred_prob_multi <- predict(fit_multi, type = "probs")
head(pred_prob_multi)

##   having fun improving    winning   relaxing
## 1  0.4898887 0.3434585 0.11961237 0.04704041
## 2  0.4909558 0.2256404 0.19476332 0.08864043
## 3  0.6343471 0.1866829 0.07281431 0.10615567
## 4  0.4655270 0.3430497 0.14620619 0.04521714
## 5  0.2953007 0.4238608 0.20265543 0.07818311
## 6  0.5219847 0.3038631 0.11767645 0.05647574

gaming_pred_multi <- gaming_clean %>%
  mutate(prediksi = predict(fit_multi, type = "class"))

conf_multi <- table(Aktual = gaming_pred_multi$whyplay, Prediksi = gaming_pred_multi$prediksi)
conf_multi

##             Prediksi
## Aktual       having fun improving winning relaxing
##   having fun       3600      1456       3        0
##   improving        2660      2003       8        0
##   winning          1019       931      12        0
##   relaxing          397       221       1        0

accuracy_multi <- sum(diag(conf_multi)) / sum(conf_multi)
accuracy_multi

## [1] 0.4560962

Interpretasi: Akurasi klasifikasi model multinomial adalah 45.61%. Nilai ini lebih tinggi dibandingkan akurasi baseline jika seluruh observasi diprediksi masuk ke kategori terbanyak, menunjukkan prediktor memberikan informasi tambahan.

1.8 Visualisasi Probabilitas Prediksi

grid_multi <- expand.grid(
  GAD_T  = seq(min(gaming_clean$GAD_T), max(gaming_clean$GAD_T), length.out = 100),
  SWL_T  = mean(gaming_clean$SWL_T),
  SPIN_T = mean(gaming_clean$SPIN_T),
  Hours  = mean(gaming_clean$Hours),
  Age    = mean(gaming_clean$Age),
  Gender = "Male"
)

grid_prob_multi <- predict(fit_multi, newdata = grid_multi, type = "probs")

grid_multi_plot <- grid_multi %>%
  bind_cols(as.data.frame(grid_prob_multi)) %>%
  pivot_longer(
    cols = c("having fun", "improving", "winning", "relaxing"),
    names_to = "whyplay", values_to = "probabilitas"
  )

ggplot(grid_multi_plot, aes(x = GAD_T, y = probabilitas, color = whyplay)) +
  geom_line(linewidth = 1.25) +
  scale_y_continuous(labels = percent_format()) +
  scale_color_manual(values = c("#3a7bd5", "#6db3f2", "#e05c5c", "#6dcfc2")) +
  labs(
    title = "Prediksi Probabilitas Motivasi Bermain berdasarkan Skor Kecemasan",
    subtitle = "Variabel lain ditahan pada nilai rata-rata; Gender = Male",
    x = "Skor GAD-7 (GAD_T)", y = "Probabilitas Prediksi", color = "Motivasi Bermain"
  ) +
  theme_adk()

Interpretasi grafik: Kurva menunjukkan bagaimana probabilitas prediksi setiap kategori motivasi bermain berubah seiring meningkatnya skor kecemasan, dengan variabel lain ditahan pada nilai rata-rata. Kategori dengan kurva yang meningkat tajam mengindikasikan motivasi yang semakin dominan pada skor kecemasan tinggi, sementara kategori dengan kurva menurun mengindikasikan motivasi yang semakin jarang dipilih.

1.9 Catatan Asumsi IIA (Independence of Irrelevant Alternatives)

Regresi logistik multinomial mengasumsikan IIA, yaitu rasio peluang antara dua kategori tidak dipengaruhi oleh ada/tidaknya kategori lain. Asumsi ini masuk akal untuk kategori motivasi bermain karena keempatnya relatif berdiri sendiri (tidak ada kategori yang merupakan subtitusi sangat dekat satu sama lain), namun tetap perlu diperhatikan secara konseptual saat menambahkan kategori motivasi baru di masa depan.

Bagian II — Regresi Logistik Ordinal

2.1 Definisi dan Tujuan Kasus

Variabel respons yang digunakan adalah GAD_Level, yaitu rekategorisasi skor GAD_T menjadi 4 tingkat keparahan kecemasan yang berurutan: Minimal < Mild < Moderate < Severe. Prediktor yang digunakan adalah kepuasan hidup (SWL_T), fobia sosial (SPIN_T), jam bermain (Hours), usia (Age), dan jenis kelamin (Gender).

\[\text{logit}\big[P(Y \le j)\big] = \theta_j - (\beta_1\,SWL\_T + \beta_2\,SPIN\_T + \beta_3\,Hours + \beta_4\,Age + \beta_5\,Gender), \quad j = \text{Minimal, Mild, Moderate}\]

2.2 Distribusi Respons Ordinal

gaming_clean %>%
  count(GAD_Level) %>%
  mutate(proporsi = round(n / sum(n), 4)) %>%
  kable(
    caption = "Tabel 3.4 Distribusi Tingkat Keparahan Kecemasan (GAD_Level)",
    col.names = c("Kategori", "Frekuensi", "Proporsi")
  ) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = FALSE
  )

Tabel 3.4 Distribusi Tingkat Keparahan Kecemasan (GAD_Level)
Kategori	Frekuensi	Proporsi
Minimal	6826	0.5545
Mild	3365	0.2733
Moderate	1466	0.1191
Severe	654	0.0531

gaming_clean %>%
  count(GAD_Level) %>%
  mutate(proporsi = n / sum(n)) %>%
  ggplot(aes(x = GAD_Level, y = proporsi, fill = GAD_Level)) +
  geom_col(width = 0.62) +
  geom_text(aes(label = percent(proporsi, accuracy = 0.1)),
            vjust = -0.4, fontface = "bold", color = "#1e3a5f") +
  scale_y_continuous(labels = percent_format(), limits = c(0, 0.5)) +
  scale_fill_manual(values = c("#6dcfc2", "#6db3f2", "#3a7bd5", "#e05c5c")) +
  labs(
    title = "Distribusi Tingkat Keparahan Kecemasan",
    subtitle = "Respons ordinal: kategori memiliki urutan alami.",
    x = NULL, y = "Proporsi"
  ) +
  theme_adk() +
  theme(legend.position = "none")

2.3 Estimasi Model Ordinal dengan R

fit_ord <- MASS::polr(
  GAD_Level ~ SWL_T + SPIN_T + Hours + Age + Gender,
  data   = gaming_clean,
  method = "logistic",
  Hess   = TRUE
)

summary(fit_ord)

## Call:
## MASS::polr(formula = GAD_Level ~ SWL_T + SPIN_T + Hours + Age + 
##     Gender, data = gaming_clean, Hess = TRUE, method = "logistic")
## 
## Coefficients:
##                  Value Std. Error t value
## SWL_T        -0.088037   0.002847 -30.921
## SPIN_T        0.052558   0.001485  35.385
## Hours         0.003559   0.001393   2.556
## Age          -0.015162   0.005914  -2.564
## GenderFemale  0.602573   0.077956   7.730
## GenderOther   0.756471   0.305909   2.473
## 
## Intercepts:
##                 Value    Std. Error t value 
## Minimal|Mild     -0.6764   0.1527    -4.4288
## Mild|Moderate     1.0130   0.1533     6.6085
## Moderate|Severe   2.5499   0.1574    16.1965
## 
## Residual Deviance: 23459.66 
## AIC: 23477.66

2.4 Ringkasan Koefisien Ordinal

ord_coef <- as.data.frame(coef(summary(fit_ord)))

result_ord <- ord_coef %>%
  rownames_to_column("parameter") %>%
  dplyr::rename(estimate = Value, std_error = `Std. Error`, t_value = `t value`) %>%
  mutate(
    p_value = 2 * (1 - pnorm(abs(t_value))),
    jenis   = ifelse(grepl("\\|", parameter), "Cutpoint", "Koefisien"),
    OR      = ifelse(jenis == "Koefisien", exp(estimate), NA_real_),
    CI_low  = ifelse(jenis == "Koefisien", exp(estimate - 1.96 * std_error), NA_real_),
    CI_high = ifelse(jenis == "Koefisien", exp(estimate + 1.96 * std_error), NA_real_)
  )

result_ord %>%
  mutate(across(c(estimate, std_error, t_value, p_value, OR, CI_low, CI_high), ~ round(.x, 4))) %>%
  kable(
    caption = "Tabel 3.5 Ringkasan Hasil Regresi Logistik Ordinal",
    col.names = c("Parameter", "Estimate", "SE", "t-value", "p-value", "Jenis", "OR", "CI 2.5%", "CI 97.5%")
  ) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = TRUE
  )

Tabel 3.5 Ringkasan Hasil Regresi Logistik Ordinal
Parameter	Estimate	SE	t-value	p-value	Jenis	OR	CI 2.5%	CI 97.5%
SWL_T	-0.0880	0.0028	-30.9212	0.0000	Koefisien	0.9157	0.9106	0.9209
SPIN_T	0.0526	0.0015	35.3853	0.0000	Koefisien	1.0540	1.0509	1.0570
Hours	0.0036	0.0014	2.5560	0.0106	Koefisien	1.0036	1.0008	1.0063
Age	-0.0152	0.0059	-2.5636	0.0104	Koefisien	0.9850	0.9736	0.9964
GenderFemale	0.6026	0.0780	7.7296	0.0000	Koefisien	1.8268	1.5680	2.1284
GenderOther	0.7565	0.3059	2.4729	0.0134	Koefisien	2.1307	1.1699	3.8808
Minimal\|Mild	-0.6764	0.1527	-4.4288	0.0000	Cutpoint	NA	NA	NA
Mild\|Moderate	1.0130	0.1533	6.6085	0.0000	Cutpoint	NA	NA	NA
Moderate\|Severe	2.5499	0.1574	16.1965	0.0000	Cutpoint	NA	NA	NA

2.5 Interpretasi Hasil Model Ordinal

Interpretasi SPIN_T: Koefisien sebesar 0.0526 (OR = 1.054) menunjukkan bahwa setiap kenaikan satu poin skor fobia sosial meningkatkan peluang berada pada kategori keparahan kecemasan yang lebih tinggi, dengan variabel lain dikontrol.

Interpretasi SWL_T: Koefisien sebesar -0.088 (OR = 0.9157) menunjukkan bahwa kepuasan hidup yang lebih tinggi berasosiasi dengan peluang yang lebih rendah untuk berada pada kategori kecemasan yang lebih parah — sesuai dugaan teoritis.

Interpretasi Hours: Koefisien sebesar 0.0036 (OR = 1.0036), p-value = 0.0106. Durasi bermain berasosiasi signifikan dengan tingkat keparahan kecemasan.

2.6 Prediksi Probabilitas dan Confusion Matrix

pred_prob_ord <- predict(fit_ord, type = "probs")
head(pred_prob_ord)

##     Minimal      Mild   Moderate      Severe
## 1 0.8039724 0.1529507 0.03348880 0.009588136
## 2 0.3991353 0.3833752 0.16108858 0.056400915
## 3 0.2836403 0.3983537 0.22686733 0.091138741
## 4 0.6446099 0.2630051 0.07096371 0.021421321
## 5 0.5224201 0.3331738 0.10938059 0.035025484
## 6 0.6193083 0.2787688 0.07809772 0.023825185

gaming_pred_ord <- gaming_clean %>%
  mutate(prediksi = predict(fit_ord, type = "class"))

conf_ord <- table(Aktual = gaming_pred_ord$GAD_Level, Prediksi = gaming_pred_ord$prediksi)
conf_ord

##           Prediksi
## Aktual     Minimal Mild Moderate Severe
##   Minimal     6055  698       70      3
##   Mild        2271  936      142     16
##   Moderate     656  587      180     43
##   Severe       180  284      120     70

accuracy_ord <- sum(diag(conf_ord)) / sum(conf_ord)
accuracy_ord

## [1] 0.5881732

Interpretasi: Akurasi klasifikasi model ordinal adalah 58.82%. Confusion matrix menunjukkan bahwa model cenderung memprediksi kategori dengan frekuensi terbanyak lebih akurat dibandingkan kategori dengan frekuensi kecil, sebuah pola umum pada data dengan distribusi kategori yang tidak seimbang.

2.7 Visualisasi Prediksi Probabilitas Ordinal

grid_ord <- expand.grid(
  SWL_T  = mean(gaming_clean$SWL_T),
  SPIN_T = seq(min(gaming_clean$SPIN_T), max(gaming_clean$SPIN_T), length.out = 100),
  Hours  = mean(gaming_clean$Hours),
  Age    = mean(gaming_clean$Age),
  Gender = "Male"
)

grid_prob_ord <- predict(fit_ord, newdata = grid_ord, type = "probs")

grid_ord_plot <- grid_ord %>%
  bind_cols(as.data.frame(grid_prob_ord)) %>%
  pivot_longer(
    cols = c("Minimal", "Mild", "Moderate", "Severe"),
    names_to = "GAD_Level", values_to = "probabilitas"
  )

ggplot(grid_ord_plot, aes(x = SPIN_T, y = probabilitas, color = GAD_Level)) +
  geom_line(linewidth = 1.25) +
  scale_y_continuous(labels = percent_format()) +
  scale_color_manual(values = c("#6dcfc2", "#6db3f2", "#3a7bd5", "#e05c5c")) +
  labs(
    title = "Prediksi Probabilitas Tingkat Kecemasan berdasarkan Skor Fobia Sosial",
    subtitle = "Variabel lain ditahan pada nilai rata-rata; Gender = Male",
    x = "Skor SPIN (SPIN_T)", y = "Probabilitas Prediksi", color = "Tingkat Kecemasan"
  ) +
  theme_adk()

Interpretasi grafik: Semakin tinggi skor fobia sosial, probabilitas berada pada kategori kecemasan Severe cenderung meningkat, sementara probabilitas kategori Minimal menurun — konsisten dengan pola cumulative logit pada model ordinal.

2.8 Pemeriksaan Asumsi Proportional Odds

# install.packages("ordinal")
library(ordinal)

fit_clm <- ordinal::clm(
  GAD_Level ~ SWL_T + SPIN_T + Hours + Age + Gender,
  data = gaming_clean,
  link = "logit"
)

nom_test <- ordinal::nominal_test(fit_clm)
nom_test

## Tests of nominal effects
## 
## formula: GAD_Level ~ SWL_T + SPIN_T + Hours + Age + Gender
##        Df logLik   AIC     LRT  Pr(>Chi)    
## <none>    -11730 23478                      
## SWL_T   2 -11718 23458 23.7990 6.794e-06 ***
## SPIN_T  2 -11730 23481  0.4762  0.788114    
## Hours   2 -11724 23471 10.6989  0.004751 ** 
## Age     2 -11728 23477  4.4300  0.109152    
## Gender  4 -11727 23480  5.7335  0.219954    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Interpretasi: Uji nominal_test membandingkan model proportional odds dengan model yang memperbolehkan efek prediktor berbeda antar batas kategori (non-proportional). Prediktor dengan p-value < 0.05 pada uji ini mengindikasikan pelanggaran asumsi proportional odds untuk variabel tersebut, sehingga interpretasi koefisiennya perlu dilakukan lebih hati-hati (atau dipertimbangkan model partial proportional odds). Prediktor dengan p-value > 0.05 dianggap memenuhi asumsi proportional odds.

Bagian III — Regresi Poisson

3.1 Definisi dan Tujuan Kasus

Variabel respons yang digunakan adalah Hours, yaitu jumlah jam bermain game per minggu — sebuah data hitung (count) nonnegatif. Prediktor yang digunakan adalah skor kecemasan (GAD_T), kepuasan hidup (SWL_T), fobia sosial (SPIN_T), usia (Age), dan jenis kelamin (Gender).

\[\ln\big[E(Hours \mid X)\big] = \beta_0 + \beta_1\,GAD\_T + \beta_2\,SWL\_T + \beta_3\,SPIN\_T + \beta_4\,Age + \beta_5\,Gender\]

Catatan: Model ini tidak menyertakan offset (exposure) karena data tidak memuat variabel periode observasi yang berbeda antar responden — seluruh responden ditanya dengan kerangka waktu yang sama (jam per minggu), berbeda dengan contoh pada materi PDF yang menggunakan offset(log(exposure)).

3.2 Distribusi Respons Hitungan

ggplot(gaming_clean, aes(x = Hours)) +
  geom_histogram(binwidth = 5, fill = "#3a7bd5", color = "white", alpha = 0.92) +
  labs(
    title = "Distribusi Jam Bermain Game per Minggu",
    subtitle = "Respons Poisson: data hitung nonnegatif.",
    x = "Jam bermain (per minggu)", y = "Frekuensi"
  ) +
  theme_adk()

3.3 Estimasi Model Poisson dengan R

fit_pois <- glm(
  Hours ~ GAD_T + SWL_T + SPIN_T + Age + Gender,
  data   = gaming_clean,
  family = poisson(link = "log")
)

summary(fit_pois)

## 
## Call:
## glm(formula = Hours ~ GAD_T + SWL_T + SPIN_T + Age + Gender, 
##     family = poisson(link = "log"), data = gaming_clean)
## 
## Coefficients:
##                Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   3.5178601  0.0158456 222.009   <2e-16 ***
## GAD_T         0.0041480  0.0004812   8.620   <2e-16 ***
## SWL_T        -0.0091438  0.0003002 -30.463   <2e-16 ***
## SPIN_T        0.0021277  0.0001639  12.980   <2e-16 ***
## Age          -0.0160653  0.0006295 -25.520   <2e-16 ***
## GenderFemale -0.1721878  0.0093263 -18.463   <2e-16 ***
## GenderOther   0.1615351  0.0294200   5.491    4e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 90913  on 12310  degrees of freedom
## Residual deviance: 87607  on 12304  degrees of freedom
## AIC: 145739
## 
## Number of Fisher Scoring iterations: 5

3.4 Ringkasan Koefisien (IRR dan Interval Kepercayaan)

pois_coef <- as.data.frame(coef(summary(fit_pois))) %>%
  rownames_to_column("parameter") %>%
  dplyr::rename(estimate = Estimate, std_error = `Std. Error`, z_value = `z value`, p_value = `Pr(>|z|)`) %>%
  mutate(
    IRR = exp(estimate),
    CI_low  = exp(estimate - 1.96 * std_error),
    CI_high = exp(estimate + 1.96 * std_error),
    perubahan_persen = 100 * (IRR - 1)
  )

pois_coef %>%
  mutate(across(c(estimate, std_error, z_value, p_value, IRR, CI_low, CI_high, perubahan_persen), ~ round(.x, 4))) %>%
  kable(
    caption = "Tabel 3.6 Ringkasan Hasil Regresi Poisson",
    col.names = c("Parameter", "Estimate", "SE", "z-value", "p-value", "IRR", "CI 2.5%", "CI 97.5%", "Perubahan (%)")
  ) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = TRUE
  )

Tabel 3.6 Ringkasan Hasil Regresi Poisson
Parameter	Estimate	SE	z-value	IRR	CI 2.5%	CI 97.5%	Perubahan (%)
(Intercept)	3.5179	0.0158	222.0093	33.7122	32.6813	34.7756	3271.2210
GAD_T	0.0041	0.0005	8.6201	1.0042	1.0032	1.0051	0.4157
SWL_T	-0.0091	0.0003	-30.4635	0.9909	0.9903	0.9915	-0.9102
SPIN_T	0.0021	0.0002	12.9805	1.0021	1.0018	1.0025	0.2130
Age	-0.0161	0.0006	-25.5198	0.9841	0.9828	0.9853	-1.5937
GenderFemale	-0.1722	0.0093	-18.4626	0.8418	0.8266	0.8574	-15.8179
GenderOther	0.1615	0.0294	5.4907	1.1753	1.1095	1.2451	17.5314

3.5 Interpretasi Koefisien

Interpretasi GAD_T: IRR = 1.0042, artinya setiap kenaikan satu poin skor kecemasan meningkatkan rata-rata jam bermain sekitar 0.42%, dengan variabel lain dikontrol.

Interpretasi SWL_T: IRR = 0.9909, artinya setiap kenaikan satu poin skor kepuasan hidup menurunkan rata-rata jam bermain sekitar 0.91%.

Interpretasi SPIN_T: IRR = 1.0021, artinya setiap kenaikan satu poin skor fobia sosial berasosiasi dengan kenaikan rata-rata jam bermain sekitar 0.21%.

3.6 Pemeriksaan Overdispersion

dispersion_pois <- sum(residuals(fit_pois, type = "pearson")^2) / df.residual(fit_pois)

tibble::tibble(
  `Dispersion Pearson` = round(dispersion_pois, 3),
  `Interpretasi ringkas` = dplyr::case_when(
    dispersion_pois < 1.5 ~ "Tidak ada indikasi overdispersion berat",
    dispersion_pois < 2.5 ~ "Ada indikasi overdispersion sedang",
    TRUE ~ "Ada indikasi overdispersion kuat"
  )
) %>%
  kable(caption = "Tabel 3.7 Indikasi Overdispersion pada Model Poisson") %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = FALSE
  )

Tabel 3.7 Indikasi Overdispersion pada Model Poisson
Dispersion Pearson	Interpretasi ringkas
7.861	Ada indikasi overdispersion kuat

Interpretasi: Rasio dispersi Pearson sebesar 7.861. Nilai ini jauh di atas 1, mengindikasikan overdispersion — variansi data hitung lebih besar daripada yang diasumsikan distribusi Poisson (rata-rata = varians). Model Negative Binomial lebih tepat digunakan sebagai alternatif.

fit_nb <- MASS::glm.nb(
  Hours ~ GAD_T + SWL_T + SPIN_T + Age + Gender,
  data = gaming_clean
)

summary(fit_nb)

## 
## Call:
## MASS::glm.nb(formula = Hours ~ GAD_T + SWL_T + SPIN_T + Age + 
##     Gender, data = gaming_clean, init.theta = 3.478783608, link = log)
## 
## Coefficients:
##                Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   3.5095688  0.0411633  85.260  < 2e-16 ***
## GAD_T         0.0038903  0.0013111   2.967   0.0030 ** 
## SWL_T        -0.0089144  0.0008004 -11.138  < 2e-16 ***
## SPIN_T        0.0021160  0.0004452   4.752 2.01e-06 ***
## Age          -0.0158259  0.0016050  -9.861  < 2e-16 ***
## GenderFemale -0.1656267  0.0236023  -7.017 2.26e-12 ***
## GenderOther   0.1844673  0.0871644   2.116   0.0343 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Negative Binomial(3.4788) family taken to be 1)
## 
##     Null deviance: 13301  on 12310  degrees of freedom
## Residual deviance: 12845  on 12304  degrees of freedom
## AIC: 94071
## 
## Number of Fisher Scoring iterations: 1
## 
## 
##               Theta:  3.4788 
##           Std. Err.:  0.0507 
## 
##  2 x log-likelihood:  -94055.1380

AIC(fit_pois, fit_nb) %>%
  tibble::rownames_to_column("model") %>%
  kable(
    caption = "Tabel 3.8 Perbandingan AIC Model Poisson vs Negative Binomial",
    col.names = c("Model", "df", "AIC")
  ) %>%
  kable_styling(
    bootstrap_options = c("striped", "hover", "condensed", "responsive"),
    full_width = FALSE
  )

Tabel 3.8 Perbandingan AIC Model Poisson vs Negative Binomial
Model	df	AIC
fit_pois	7	145739.05
fit_nb	8	94071.14

Interpretasi: Model Negative Binomial memiliki AIC yang lebih rendah dibandingkan model Poisson, sehingga lebih disarankan untuk data ini karena mampu mengakomodasi overdispersion pada variabel Hours.

3.7 Visualisasi Prediksi Rate Poisson

grid_pois <- expand.grid(
  GAD_T  = seq(min(gaming_clean$GAD_T), max(gaming_clean$GAD_T), length.out = 100),
  SWL_T  = mean(gaming_clean$SWL_T),
  SPIN_T = mean(gaming_clean$SPIN_T),
  Age    = mean(gaming_clean$Age),
  Gender = "Male"
)

pred_pois <- predict(fit_pois, newdata = grid_pois, type = "link", se.fit = TRUE)

grid_pois_plot <- grid_pois %>%
  mutate(
    fit_link = pred_pois$fit,
    se_link  = pred_pois$se.fit,
    rate     = exp(fit_link),
    lower    = exp(fit_link - 1.96 * se_link),
    upper    = exp(fit_link + 1.96 * se_link)
  )

ggplot(grid_pois_plot, aes(x = GAD_T, y = rate)) +
  geom_ribbon(aes(ymin = lower, ymax = upper), fill = "#bdd7ff", alpha = 0.6) +
  geom_line(color = "#3a7bd5", linewidth = 1.25) +
  labs(
    title = "Prediksi Rata-rata Jam Bermain berdasarkan Skor Kecemasan",
    subtitle = "Variabel lain ditahan pada nilai rata-rata; Gender = Male; pita biru = interval kepercayaan 95%",
    x = "Skor GAD-7 (GAD_T)", y = "Prediksi Jam Bermain per Minggu"
  ) +
  theme_adk()

Interpretasi grafik: Garis menunjukkan rata-rata prediksi jam bermain per minggu pada berbagai tingkat skor kecemasan, dengan pita biru sebagai interval kepercayaan 95%. Tren garis yang naik mengonfirmasi asosiasi positif antara skor kecemasan dan durasi bermain game.

🎯 Kesimpulan Tugas 3

Regresi Logistik Multinomial — Model dengan respons whyplay (4 kategori motivasi bermain) menghasilkan uji rasio kemungkinan G² = 662.3795 (df = 21, p-value = 0), sehingga prediktor secara keseluruhan signifikan menjelaskan variasi motivasi bermain game. Akurasi klasifikasi model sebesar 45.61%.
Regresi Logistik Ordinal — Model dengan respons GAD_Level (4 tingkat keparahan kecemasan) menunjukkan bahwa skor fobia sosial (SPIN_T) memiliki OR = 1.054, mengindikasikan bahwa fobia sosial yang lebih tinggi berasosiasi dengan tingkat kecemasan yang lebih parah. Akurasi klasifikasi model sebesar 58.82%.
Regresi Poisson — Model dengan respons Hours (jam bermain per minggu) menghasilkan rasio dispersi Pearson sebesar 7.861, yang mengindikasikan overdispersion sehingga model Negative Binomial menjadi alternatif yang lebih sesuai dibandingkan Poisson murni.
Kesimpulan umum: Ketiga model GLM ini menunjukkan bahwa profil psikologis pemain — kecemasan (GAD-7), kepuasan hidup (SWL), dan fobia sosial (SPIN) — berasosiasi dengan motivasi bermain, tingkat keparahan kecemasan, dan intensitas bermain game secara konsisten, sesuai dengan jenis skala pengukuran masing-masing variabel respons (nominal, ordinal, dan hitungan).

Referensi

Agresti, A. (2013). Categorical Data Analysis (3rd ed.). Wiley.
Mindra Jaya, I. G. N. (2024). Analisis Data Kategori. RPubs. https://rpubs.com/mindra/DataKategori
R Core Team (2024). R: A Language and Environment for Statistical Computing. R Foundation.
Wickham, H. (2016). ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag.

Laporan ini dibuat menggunakan R Markdown · Dipublikasikan ke RPubs · Hanifah Putri Hariel · S1 Statistika FMIPA Unpad

📊 Analisis Data Kategori

Hanifah Putri Hariel

2026-05-31

Mata Kuliah: Analisis Data Kategori

Tugas 1 — Konsep Dasar Analisis Data Kategori

1. Definisi Analisis Data Kategori

2. Tabel Kontingensi

2.1 Definisi

2.2 Struktur Tabel Kontingensi 2×2

2.3 Konsep Distribusi

3. Ukuran Asosiasi

3.1 Odds

3.2 Odds Ratio (OR)

3.3 Relative Risk (RR)

4. Contoh Perhitungan Manual

4.1 Tabel Kontingensi

4.2 Menghitung Peluang Bersyarat

4.3 Menghitung Odds

4.4 Menghitung Odds Ratio

5. Analisis Menggunakan R

6. Interpretasi Hasil

Tugas 6 — Inferensi Tabel Kontingensi Dua Arah

Pendahuluan

Kasus 1 — Tabel Kontingensi 2×2: Merokok dan Kanker Paru

1.1 Data dan Penyusunan Tabel

1.2 Estimasi Titik Proporsi

1.3 Interval Kepercayaan 95%

Proporsi Masing-Masing Kelompok

Risk Difference (RD)

Relative Risk (RR)

Odds Ratio (OR)

Ringkasan Ukuran Asosiasi

1.4 Uji Dua Proporsi

1.5 Uji Chi-Square Independensi

1.6 Uji Likelihood Ratio (G²)

1.7 Fisher Exact Test

1.8 Perbandingan Keempat Uji

1.9 Visualisasi Kasus 1

🎯 Kesimpulan Kasus 1

Kasus 2 — Tabel Kontingensi 2×3: Gender dan Identifikasi Partai Politik

2.1 Data dan Penyusunan Tabel

2.2 Frekuensi Harapan

2.3 Uji Chi-Square Independensi

2.4 Residual Pearson dan Standardized Residual

2.5 Partisi Chi-Square

Partisi 1 — Democrat vs Republican

Partisi 2 — (Democrat + Republican) vs Independent

Perbandingan Partisi dengan Uji Keseluruhan

2.6 Visualisasi Kasus 2

2.7 Identifikasi Kategori yang Paling Berkontribusi

🎯 Kesimpulan Kasus 2

Tugas 3 — Praktikum Regresi Logistik Multinomial, Ordinal, dan Regresi Poisson

Pendahuluan

Definisi Singkat

Persiapan: Package dan Data

Instalasi dan Pemanggilan Package

Input Data

Variabel yang Digunakan

Pembersihan dan Persiapan Data

Bagian I — Regresi Logistik Multinomial

1.1 Definisi dan Tujuan Kasus

1.2 Distribusi Kategori Respons

1.3 Estimasi Model Multinomial dengan R

1.4 Ringkasan Koefisien, Standard Error, z-value, p-value, dan RRR

1.5 Interpretasi Koefisien Multinomial

1.6 Uji Signifikansi Model Secara Keseluruhan (Likelihood Ratio Test)

1.7 Prediksi Probabilitas dan Confusion Matrix

1.8 Visualisasi Probabilitas Prediksi

1.9 Catatan Asumsi IIA (Independence of Irrelevant Alternatives)

Bagian II — Regresi Logistik Ordinal

2.1 Definisi dan Tujuan Kasus

2.2 Distribusi Respons Ordinal

2.3 Estimasi Model Ordinal dengan R

2.4 Ringkasan Koefisien Ordinal

2.5 Interpretasi Hasil Model Ordinal

2.6 Prediksi Probabilitas dan Confusion Matrix

2.7 Visualisasi Prediksi Probabilitas Ordinal

2.8 Pemeriksaan Asumsi Proportional Odds

Bagian III — Regresi Poisson

3.1 Definisi dan Tujuan Kasus