Pendahuluan

Data kategori, juga dikenal sebagai data kualitatif, adalah jenis data yang mengelompokkan observasi ke dalam satu atau lebih kategori berdasarkan karakteristik atau atribut, bukan nilai numerik. Menurut Statgraphics, data kategori mengklasifikasikan observasi ke dalam kategori, seperti “setuju” atau “tidak setuju” dalam survei.

Analisis data kategori adalah proses menggunakan metode statistik untuk menganalisis, menginterpretasikan, dan memahami data yang dikelompokkan ke dalam kategori. Menurut StudySmarter, analisis ini melibatkan pemeriksaan, interpretasi, dan presentasi data yang dikelompokkan secara kualitatif ke dalam kelompok ordinal atau nominal.

Tujuan Analisis Data Kategori

Mengelompokkan Data

Tujuan pertama adalah mengorganisasi data ke dalam kategori yang bermakna berdasarkan karakteristik tertentu. Ini melibatkan pengelompokan observasi ke dalam kelompok seperti gender (laki-laki, perempuan), status perkawinan (lajang, menikah, bercerai), atau preferensi produk (elektronik, pakaian, makanan). Langkah ini sangat penting untuk menyusun data agar memudahkan analisis lebih lanjut. Misalnya, dalam penelitian pasar, mengelompokkan pelanggan berdasarkan rentang usia (misalnya, 18-25, 26-35, dll.) membantu dalam segmentasi pasar untuk kampanye yang ditargetkan.

Memahami Pola dan Hubungan

Tujuan kunci lainnya adalah mengidentifikasi pola, tren, dan korelasi dalam dan antar variabel kategori. Ini melibatkan analisis bagaimana kategori saling berhubungan, seperti apakah ada hubungan antara gender dan preferensi produk tertentu. Metode seperti tabel frekuensi, grafik batang, dan grafik lingkaran sering digunakan untuk memvisualisasikan distribusi, sementara uji statistik seperti uji chi-kuadrat atau uji Fisher’s exact digunakan untuk menilai hubungan.

Dalam ilmu sosial, peneliti mungkin menggunakan analisis data kategori untuk menentukan apakah ada korelasi antara tingkat pendidikan (misalnya, SMA, kuliah, pascasarjana) dan preferensi partai politik.

Mendukung Pengambilan Keputusan

Hasil analisis data kategori sering digunakan untuk mendukung proses pengambilan keputusan di berbagai sektor. Dalam bisnis, misalnya, perusahaan mungkin menganalisis data kategori untuk mengelompokkan pelanggan dan menyesuaikan strategi pemasaran, seperti menargetkan iklan ke kelompok demografis tertentu berdasarkan preferensi pembelian. Dalam kebijakan publik, pemerintah mungkin menggunakan analisis ini untuk mengevaluasi dampak kebijakan pada segmen populasi yang berbeda, seperti menilai efektivitas program kesejahteraan sosial untuk kelompok berpenghasilan rendah versus tinggi.

Rantai ritel mungkin menggunakan data kategori untuk memutuskan produk mana yang akan disediakan di wilayah berbeda berdasarkan preferensi pelanggan (misalnya, perkotaan vs. pedesaan).

Mengembangkan Model Prediktif

Analisis data kategori juga bertujuan mengembangkan model prediktif yang meramalkan hasil berdasarkan variabel kategori. Teknik seperti regresi logistik sering digunakan, terutama untuk hasil biner (misalnya, ya/tidak, sukses/gagal). Misalnya, dalam kesehatan, sebuah model mungkin memprediksi kemungkinan pasien mengalami penyakit berdasarkan faktor kategori seperti status merokok (ya/tidak) atau riwayat keluarga (ya/tidak). Tujuan ini sangat penting untuk perencanaan strategis dan peramalan di bidang seperti bisnis dan kedokteran.

Dalam pemasaran, model prediktif mungkin meramalkan churn pelanggan berdasarkan variabel kategori seperti jenis langganan (bulanan, tahunan) dan kategori umpan balik pelanggan (puas, tidak puas).

Definisi dan Ruang Lingkup Analisis Data Kategori

Definisi Analisis Data Kategori

Analisis data kategori didefinisikan sebagai serangkaian teknik statistik yang digunakan untuk menganalisis, menginterpretasikan, dan menyajikan data yang telah dikelompokkan ke dalam kategori berdasarkan karakteristik kualitatif. Menurut Statgraphics, analisis ini berfokus pada data yang tidak dapat diukur secara kuantitatif tetapi dapat diklasifikasikan, seperti data nominal (kategori tanpa urutan alami) dan data ordinal (kategori dengan urutan).

  • Data Nominal: Merujuk pada kategori yang tidak memiliki urutan atau hierarki alami. Contohnya adalah jenis kelamin (laki-laki, perempuan), warna (merah, biru, hijau), atau status perkawinan (lajang, menikah, bercerai). Data ini hanya mengidentifikasi perbedaan tanpa menunjukkan tingkat atau urutan.

  • Data Ordinal: Merujuk pada kategori yang memiliki urutan atau tingkat, tetapi jarak antar kategori tidak selalu sama. Contohnya adalah tingkat pendidikan (SD, SMP, SMA, S1), tingkat kepuasan (sangat puas, puas, netral, tidak puas), atau skala Likert (1 hingga 5). Meskipun ada urutan, perbedaan kuantitatif antar kategori tidak dapat dihitung dengan tepat.

Ruang Lingkup Analisis Data Kategori

Ruang lingkup analisis data kategori mencakup berbagai aspek klasifikasi data dan penerapannya dalam konteks praktis. Berikut adalah penjelasan rinci berdasarkan subtopik yang relevan:

  1. Nominal vs Ordinal

    • Nominal: Seperti disebutkan, data nominal tidak memiliki urutan. Analisisnya biasanya terbatas pada perhitungan frekuensi atau proporsi dalam setiap kategori. Misalnya, dalam survei, kita mungkin menghitung persentase responden yang memilih “ya” atau “tidak”

    • Ordinal: Data ordinal memungkinkan analisis yang lebih mendalam karena adanya urutan. Misalnya, kita dapat menilai tren, seperti apakah tingkat kepuasan meningkat seiring dengan peningkatan pelayanan. Namun, karena jarak antar kategori tidak sama, penggunaan statistik parametrik seperti rata-rata harus dihindari; alih-alih, median atau modus lebih tepat.

    • Perbedaan ini penting karena memengaruhi metode analisis yang dipilih. Formpl.us menjelaskan bahwa pemahaman ini membantu dalam memilih alat statistik yang sesuai, seperti uji Mann-Whitney untuk data ordinal dibandingkan uji chi-kuadrat untuk data nominal.

  2. Data Biner vs Multikategori

  • Data Biner: Merujuk pada data yang hanya memiliki dua kategori, seperti ya/tidak, benar/salah, atau 0/1. Analisis data biner sering digunakan dalam studi kasus-kontrol atau eksperimen sederhana, seperti mengevaluasi apakah pasien sembuh (ya/tidak) setelah perawatan tertentu. Teknik seperti regresi logistik sering diterapkan untuk memodelkan hubungan antara variabel biner dan prediktor.

  • Data Multikategori: Merujuk pada data dengan lebih dari dua kategori, seperti jenis darah (A, B, AB, O), hari dalam seminggu (Senin, Selasa, dll.), atau preferensi makanan (sayuran, daging, ikan). Analisis ini lebih kompleks dan mungkin memerlukan pendekatan seperti analisis korespondensi atau uji chi-kuadrat dengan derajat bebas yang lebih tinggi.

  1. Penerapan di Berbagai Bidang

Ruang lingkup analisis data kategori juga mencakup penerapannya di berbagai disiplin ilmu, seperti:

  • Ilmu Sosial dan Psikologi: Menganalisis data seperti preferensi politik berdasarkan kategori demografi.

  • Kesehatan dan Kedokteran: Mengevaluasi hubungan antara kategori seperti status merokok dan risiko penyakit.

  • Pemasaran dan Bisnis: Mengelompokkan pelanggan berdasarkan preferensi untuk strategi pemasaran.

Metode dalam Analisis Data Kategori

Tabel Kontingensi & Uji Chi-Square

Tabel kontingensi adalah alat penting dalam analisis data kategorikal untuk merangkum hubungan antara dua atau lebih variabel kategorikal. Uji chi-square, di sisi lain, digunakan untuk menguji hubungan antar variabel dalam tabel kontingensi, khususnya untuk mengevaluasi apakah variabel-variabel tersebut independen atau memiliki asosiasi.

Tabel Kontingensi

Tabel kontingensi merangkum data kategorikal dalam bentuk frekuensi pengamatan untuk kombinasi kategori dari dua atau lebih variabel. Misalnya, tabel 2x2 menunjukkan hubungan antara dua variabel biner (masing-masing memiliki dua kategori). Struktur probabilitas untuk tabel kontingensi melibatkan :

Probabilitas Bersama (Joint Probabilities)

Probabilitas bahwa observasi berada pada kombinasi tertentu dari kategori dua variabel, misalnya \[ P(X = i, Y = j) \]

Probabilitas Marjinal (Marginal Probabilities)

Total probabilitas untuk satu variabel, dihitung dengan menjumlahkan probabilitas bersama di seluruh kategori variabel lain, misalnya \[ P(X = i) = \sum_j P(X = i, Y = j) \]

Probabilitas Bersyarat (Conditional Probabilities)

Probabilitas satu variabel diberikan nilai tertentu dari variabel lain, misalnya \[ P(Y = j | X = i) \]

Membandingkan Proporsi dalam Tabel 2x2

Dalam analisis tabel kontingensi 2x2, metode umum melibatkan perbandingan proporsi antara dua kelompok. Metode ini berguna untuk mengevaluasi hubungan antara dua variabel kategorikal.

Perbedaan Proporsi (Difference of Proportions)

Perbedaan proporsi mengukur selisih antara dua proporsi, misalnya: \[ \hat{\pi}_1 - \hat{\pi}_2 \] di mana:

  • \(\hat{\pi}_1\) adalah proporsi kejadian dalam kelompok pertama,

  • \(\hat{\pi}_2\) adalah proporsi kejadian dalam kelompok kedua.

Standar error (SE) dari perbedaan proporsi dihitung dengan rumus: \[ SE = \sqrt{\frac{\hat{\pi}_1(1-\hat{\pi}_1)}{n_1} + \frac{\hat{\pi}_2(1-\hat{\pi}_2)}{n_2}} \] di mana \(n_1\) dan \(n_2\) adalah ukuran sampel untuk masing-masing kelompok.

Contoh (Aspirin dan Serangan Jantung)

Misalkan kita memiliki data yang menunjukkan bahwa proporsi serangan jantung lebih rendah pada kelompok yang mengonsumsi aspirin dibandingkan dengan kelompok plasebo. Perbedaan proporsi dapat dihitung untuk mengevaluasi efek aspirin terhadap pengurangan serangan jantung.

Risiko Relatif (Relative Risk)

Risiko relatif mengukur rasio proporsi kejadian antara dua kelompok, yaitu: \[ \frac{\hat{\pi}_1}{\hat{\pi}_2} \] di mana:

  • \(\hat{\pi}_1\) adalah proporsi kejadian dalam kelompok pertama

  • \(\hat{\pi}_2\) adalah proporsi kejadian dalam kelompok kedua.

Risiko relatif memberikan gambaran tentang seberapa besar kemungkinan suatu kejadian terjadi di satu kelompok dibandingkan dengan kelompok lain

Rasio Odds

Rasio odds adalah ukuran asosiasi penting dalam tabel 2x2, didefinisikan sebagai rasio odds kejadian pada satu kelompok terhadap kelompok lain. Untuk tabel 2x2 dengan sel \((n_{11}, n_{12}, n_{21}, n_{22})\):

Odds Ratio: \[ \theta = \frac{n_{11}n_{22}}{n_{12}n_{21}} \]

Sifat Rasio Odds

Rasio odds bersifat simetris (tidak berubah jika baris dan kolom ditukar) dan tidak bergantung pada urutan kategori.

Contoh

Dalam studi aspirin dan serangan jantung, rasio odds menunjukkan bahwa odds serangan jantung lebih rendah pada kelompok aspirin.

Inferensi untuk Rasio Odds

Interval kepercayaan untuk rasio odds dihitung menggunakan log odds ratio, dengan standar error: \[ SE = \sqrt{\frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}}} \]

Studi Kasus-Kontrol

Rasio odds sangat relevan dalam studi kasus-kontrol karena dapat diestimasi meskipun desainnya retrospektif.

Uji Chi-Square

Uji Chi-Square dalam Tabel Kontingensi

Uji chi-square digunakan untuk menguji hipotesis nol bahwa dua variabel dalam tabel kontingensi adalah independen. Ada dua jenis statistik chi-square yang umum digunakan:

Statistik Pearson

Statistik Pearson dihitung sebagai: \[ X^2 = \sum \frac{(O_{ij} - E_{ij})^2}{E_{ij}} \] di mana: - \(O_{ij}\) adalah frekuensi observasi, - \(E_{ij} = \frac{n_{i+} n_{+j}}{n}\) adalah frekuensi yang diharapkan di bawah hipotesis independensi.

Statistik ini mengikuti distribusi chi-square dengan derajat kebebasan: \[ df = (r-1)(c-1) \] di mana \(r\) dan \(c\) adalah jumlah baris dan kolom pada tabel kontingensi.

Statistik Rasio Likelihud (Likelihood-Ratio Statistic)

Statistik rasio likelihood dihitung sebagai: \[ G^2 = 2 \sum O_{ij} \log \left( \frac{O_{ij}}{E_{ij}} \right) \] Statistik ini juga mengikuti distribusi chi-square dengan derajat kebebasan yang sama dengan statistik Pearson.

Contoh

Uji chi-square digunakan untuk menguji apakah afiliasi politik independen dari jenis kelamin. Hasil menunjukkan asosiasi yang signifikan, dengan nilai \(X^2\) dan \(G^2\) besar, yang menolak hipotesis nol.

Residu

Residus standar dihitung sebagai: \[ \frac{O_{ij} - E_{ij}}{\sqrt{E_{ij}}} \] Residus standar terstandarisasi digunakan untuk mengidentifikasi sel mana yang menyumbang paling besar terhadap asosiasi.

Partisi Chi-Square

Statistik chi-square dapat dipartisi untuk menganalisis komponen spesifik dari asosiasi. Misalnya, untuk tabel \(I \times J\), kontribusi dari sub-tabel tertentu dapat dipisahkan untuk analisis lebih lanjut.

Komentar tentang Uji Chi-Square

Uji chi-square sensitif terhadap ukuran sampel. Untuk sampel besar, bahkan asosiasi kecil dapat signifikan secara statistik, tetapi mungkin tidak signifikan secara praktis. Untuk sampel kecil, pendekatan eksak (seperti uji Fisher) lebih disarankan.

Uji Indepedensi Untuk Data Ordinal

Untuk variabel ordinal, metode khusus dapat memanfaatkan urutan kategori untuk meningkatkan daya uji:

Alternatif Tren Linear

Uji tren linear mengasumsikan hubungan linier antara skor ordinal dan proporsi. Skor (misalnya, 1, 2, 3, dst.) diberikan pada kategori, dan uji dilakukan untuk mendeteksi tren dalam proporsi.

Contoh (Penggunaan Alkohol dan Malformasi Bayi)

Uji tren linear menunjukkan hubungan antara tingkat konsumsi alkohol dan risiko malformasi bayi, dengan daya uji lebih besar dibandingkan uji chi-square nominal.

Pilihan Skor

Pemilihan skor (misalnya, equidistant atau midranks) memengaruhi hasil uji. Skor yang sesuai dengan data meningkatkan sensitivitas uji.

Uji Tren untuk Tabel \(I \times 2\) atau \(2 \times J\)

Uji ini fokus pada tren linier dalam proporsi di salah satu dimensi tabel.

Regresi Logistik

Regresi logistik adalah model statistik yang digunakan untuk menganalisis hubungan antara variabel respons kategorikal (biasanya biner, seperti “ya/tidak” atau “sukses/gagal”) dan satu atau lebih variabel prediktor (penjelas), yang dapat berupa kategorikal atau kontinu. Regresi logistik merupakan bagian dari kelas Generalized Linear Models (GLM) dan sangat penting dalam analisis data kategorikal, terutama untuk memodelkan probabilitas kejadian tertentu.

GLM terdiri dari tiga komponen utama:

  1. Komponen Acak: Untuk regresi logistik, variabel respons diasumsikan mengikuti distribusi binomial, dengan probabilitas sukses \(\pi\).

2. Komponen Sistematis: Kombinasi linear dari prediktor, yaitu: \[ \eta = \beta_0 + \beta_1x_1 + \cdots + \beta_px_p \] di mana :

  • \(x_i\) adalah prediktor dan \(\beta_i\) adalah koefisien.
3\. **Fungsi Link**: Untuk regresi logistik, fungsi link adalah logit, didefinisikan sebagai: $$
   \text{logit}(\pi) = \log\left(\frac{\pi}{1-\pi}\right) = \beta_0 + \beta_1x_1 + \cdots + \beta_px_p
   $$ Fungsi logit mengubah probabilitas $\pi$ (yang berkisar antara 0 dan 1) menjadi skala tak terbatas, memungkinkan pemodelan linier.

Model Regresi Logistik

Model ini memodelkan probabilitas sukses \(\pi\) sebagai: \[ \pi(x) = \frac{\exp(\beta_0 + \beta_1x_1 + \cdots + \beta_px_p)}{1 + \exp(\beta_0 + \beta_1x_1 + \cdots + \beta_px_p)} \]Alternatifnya, model probit menggunakan fungsi distribusi kumulatif normal standar sebagai link, tetapi regresi logistik lebih umum digunakan karena interpretasi rasio odds yang lebih mudah.

Interpretasi Model Regresi Logistik

Interpretasi parameter dalam regresi logistik sangat penting untuk memahami efek prediktor terhadap variabel respons: - Interpretasi Rasio Odds: Koefisien \(\beta_i\) dalam model logistik menunjukkan perubahan log odds untuk setiap unit kenaikan pada prediktor \(x_i\). Rasio odds dihitung sebagai: \[ \text{Odds Ratio} = \exp(\beta_i) \] Misalnya, jika \(\beta_1 = 0.693\), maka \(\exp(0.693) \approx 2\), yang berarti odds kejadian meningkat dua kali lipat untuk setiap kenaikan satu unit pada \(x_1\).

- Aproksimasi Linear : Untuk perubahan kecil pada \(x_i\), efek \(\beta_i\) dapat diaproksimasi sebagai perubahan probabilitas: \[ \frac{\partial \pi}{\partial x_i} \approx \beta_i \pi(1-\pi) \] Ini berguna untuk memahami tingkat perubahan probabilitas di sekitar nilai tertentu.

Inferensi untuk Regresi Logistik

Inferensi statistik dalam regresi logistik melibatkan estimasi parameter, pengujian signifikansi, dan pembentukan interval kepercayaan: - Estimasi Parameter: Parameter \(\beta\) diestimasi menggunakan Maximum Likelihood Estimation (MLE). Data dapat berupa kelompok (grouped) atau tidak dikelompokkan (ungrouped), tetapi hasilnya konsisten. - Interval Kepercayaan: Interval kepercayaan untuk \(\beta_i\) dihitung menggunakan pendekatan Wald: \[ \hat{\beta}_i \pm z_{\alpha/2} \cdot SE(\hat{\beta}_i) \]Untuk rasio odds, interval kepercayaan dihitung dengan mengambil eksponen dari interval untuk \(\beta_i\). - Pengujian Signifikansi: Uji Wald digunakan untuk menguji hipotesis nol \(H_0: \beta_i = 0\), dengan statistik: \[ z = \frac{\hat{\beta}_i}{SE(\hat{\beta}_i)} \]Alternatifnya, uji rasio likelihud (likelihood-ratio test) membandingkan model dengan dan tanpa prediktor tertentu.

Regresi Logistik dengan Prediktor Kategorikal

Untuk prediktor kategorikal, variabel indikator (dummy variables) digunakan untuk merepresentasikan kategori: - Variabel Indikator: Untuk prediktor dengan \(k\) kategori, \(k-1\) variabel indikator dibuat, dengan satu kategori sebagai referensi. Koefisien menunjukkan efek kategori tersebut relatif terhadap kategori referensi.

Regresi Logistik Berganda

Regresi logistik berganda melibatkan lebih dari satu prediktor: - Contoh: Model dengan dua prediktor (warna dan lebar karapas) menunjukkan efek masing-masing prediktor pada probabilitas keberadaan satelit. - Perbandingan Model: Uji rasio likelihud digunakan untuk memeriksa apakah prediktor tertentu diperlukan dalam model. - Interaksi: Termasuk istilah interaksi (misalnya, \(x_1 \cdot x_2\)) memungkinkan efek prediktor bervariasi berdasarkan nilai prediktor lain.

Membangun dan Menerapkan Model Regresi Logistik

Bab ini membahas strategi praktis dalam membangun dan memvalidasi model regresi logistik:

- Strategi Pemilihan Model: Jumlah prediktor harus dibatasi berdasarkan ukuran sampel untuk menghindari overfitting.

- Kriteria AIC: Akaike Information Criterion digunakan untuk membandingkan model, menyeimbangkan kecocokan dan kompleksitas.

- Tabel Klasifikasi dan Kurva ROC: Digunakan untuk mengevaluasi kekuatan prediktif model. Kurva ROC mengukur trade-off antara sensitivitas dan spesifisitas.

- Pemeriksaan Model: Uji Rasio Likelihud digunakan untuk menguji kecocokan model, dan deviansi digunakan untuk menguji goodness-of-fit.

Analisis Korespondensi

Analisis Korespondensi adalah metode eksplorasi yang kuat untuk memvisualisasikan hubungan antara variabel kategorikal dalam tabel kontingensi. Dengan memplot kategori baris dan kolom dalam ruang dimensi rendah, CA mengungkapkan pola asosiasi melalui kedekatan kategori dalam plot. Total inersia, yang terkait dengan statistik chi-square, didekomposisi untuk menghasilkan sumbu-sumbu yang menjelaskan variasi asosiasi. CA standar cocok untuk variabel nominal, sementara pendekatan ordinal digunakan untuk kategori terurut. Analisis Korespondensi Berganda memperluas metode ini untuk lebih dari dua variabel. Buku ini menekankan visualisasi grafis, fleksibilitas metode, dan penggunaan perangkat lunak untuk implementasi praktis, dengan contoh seperti hubungan pekerjaan dan kepuasan kerja untuk mengilustrasikan penerapan.

Analisis korespondensi bertujuan untuk merangkum hubungan dalam tabel kontingensi \(I \times J\) (dengan \(I\) baris dan \(J\) kolom) melalui visualisasi grafis yang menunjukkan kedekatan antar kategori. Metode ini berfokus pada profil baris (proporsi relatif dalam setiap baris) dan profil kolom (proporsi relatif dalam setiap kolom), serta bagaimana profil ini menyimpang dari independensi.

Tabel Kontingensi dan Profil

Untuk tabel kontingensi dengan frekuensi \(n_{ij}\), total frekuensi \(n\), dan probabilitas \(p_{ij} = \frac{n_{ij}}{n}\):

  • Profil baris: \(\frac{p_{ij}}{p_{i+}}\), di mana \(p_{i+} = \sum_j p_{ij}\) adalah total marjinal baris.
  • Profil kolom: \(\frac{p_{ij}}{p_{+j}}\), di mana \(p_{+j} = \sum_i p_{ij}\) adalah total marjinal kolom.

Jika baris dan kolom independen, profil baris akan sama untuk semua baris, dan profil kolom akan sama untuk semua kolom.

Tujuan CA

Tujuan utama analisis korespondensi adalah menemukan skor (koordinat) untuk kategori baris dan kolom dalam ruang dimensi rendah yang memaksimalkan asosiasi antar kategori, sehingga kategori yang serupa ditempatkan berdekatan dalam plot.

Dekomposisi Statistik Chi-Square

Analisis korespondensi terkait erat dengan statistik chi-square untuk menguji independensi dalam tabel kontingensi. Statistik chi-square Pearson didekomposisi untuk menghasilkan koordinat dalam CA:

  • Statistik Chi-Square: \[ X^2 = n \sum_i \sum_j \frac{(p_{ij} - p_{i+}p_{+j})^2}{p_{i+}p_{+j}} \] Statistik ini mengukur penyimpangan dari independensi. Dalam CA, \(X^2/n\) dikenal sebagai total inersia, yang mengukur total variasi dalam tabel yang dijelaskan oleh asosiasi.

  • Dekomposisi Inersia: Total inersia didekomposisi menjadi komponen-komponen yang terkait dengan sumbu-sumbu (dimensi) dalam plot CA. Setiap sumbu memiliki nilai singular (singular value) yang menunjukkan kontribusi asosiasi.

  • Jumlah Dimensi: Untuk tabel \(I \times J\), jumlah maksimum dimensi adalah \(\min(I-1, J-1)\). Biasanya, dua sumbu pertama (dimensi 1 dan 2) digunakan untuk visualisasi, karena sumbu ini menjelaskan sebagian besar inersia.

Menentuan Skor dan Plot Korespondensi

  • Skor Baris dan Kolom: CA menghasilkan skor (koordinat) untuk setiap kategori baris dan kolom menggunakan singular value decomposition (SVD) dari matriks residu terstandarisasi: \[ r_{ij} = \frac{p_{ij} - \sqrt{p_{i+}p_{+j}}}{\sqrt{p_{i+}p_{+j}}} \] Skor ini memastikan jarak antar kategori merefleksikan asosiasi. Skor baris dan kolom diplot dalam ruang dua dimensi, di mana:

    • Kategori dengan profil serupa ditempatkan berdekatan.
    • Jarak antara kategori baris dengan kategori kolom menunjukkan kekuatan asosiasi.
  • Standarisasi Skor: Skor dapat distandarisasi sehingga varians rata-rata baris dan kolom sesuai dengan kebutuhan interpretasi. Standarisasi simetris (biplot) sering digunakan untuk visualisasi.

  • Interpretasi Plot:

    • Kategori baris yang berdekatan memiliki profil kolom yang serupa, dan sebaliknya.
    • Kategori baris dan kolom yang berdekatan menunjukkan asosiasi positif yang kuat.
    • Jarak dari pusat (0,0) menunjukkan seberapa jauh kategori menyimpang dari independensi.

Asosiasi dalam Tabel dengan Skala Ordinal

Untuk tabel dengan kategori ordinal, CA dapat dimodifikasi untuk memanfaatkan urutan kategori:

  • Skor Ordinal: Skor kategori dapat dibatasi untuk mencerminkan urutan alami, menghasilkan analisis korespondensi kanonik yang lebih terfokus pada tren linier.

Analisis Korespondensi Berganda (Multiple Correspondence Analysis - MCA)

MCA adalah ekstensi dari CA untuk menganalisis hubungan antara lebih dari dua variabel kategorikal, biasanya direpresentasikan dalam tabel indikator (matriks biner yang menunjukkan keanggotaan kategori):

  • Tabel Indikator: Setiap observasi direpresentasikan sebagai vektor biner yang menunjukkan kategori untuk setiap variabel. MCA menganalisis tabel indikator ini untuk menemukan hubungan antar semua kategori.

  • Interpretasi MCA: Plot MCA menunjukkan kedekatan antar kategori dari berbagai variabel, membantu mengidentifikasi pola kompleks.

Komentar dan Aplikasi Praktis

  • Keunggulan CA:
    • Tidak memerlukan asumsi distribusi, sehingga fleksibel untuk berbagai jenis data kategorikal.
    • Visualisasi grafis memudahkan interpretasi pola asosiasi.
    • Dapat menangani tabel besar dengan banyak kategori.
  • Keterbatasan:
    • Interpretasi plot CA bersifat kualitatif dan dapat subjektif.
    • Tidak memberikan uji statistik formal untuk hipotesis spesifik.
    • Untuk tabel dengan frekuensi rendah atau sel kosong, hasil mungkin tidak stabil.
  • Aplikasi:
    • Ilmu Sosial: Menganalisis hubungan antara status sosial, pendidikan, dan preferensi politik.
    • Pemasaran: Mengidentifikasi pola pembelian berdasarkan demografi pelanggan.
    • Ekologi: Menghubungkan spesies dengan karakteristik lingkungan.

Decision Tree & Random Forest

Decision Tree

Decision Tree adalah algoritma pembelajaran mesin yang digunakan untuk tugas klasifikasi dan regresi. Dalam konteks analisis data kategorikal, Decision Tree sering digunakan untuk memprediksi variabel respons kategorikal (misalnya, “ya/tidak” atau kelas tertentu) berdasarkan variabel prediktor, yang bisa berupa kategorikal atau kontinu. Decision Tree bekerja dengan membagi data menjadi subset berdasarkan aturan keputusan yang dibentuk dari fitur prediktor, membentuk struktur seperti pohon.

Konsep Dasar Decision Tree

  • Struktur Pohon:
    • Node Akar (Root Node): Titik awal pohon, mewakili seluruh dataset.
    • Node Internal: Titik keputusan berdasarkan fitur tertentu (misalnya, “Apakah usia > 30?”).
    • Node Daun (Leaf Node): Titik akhir yang mewakili kelas atau nilai prediksi (misalnya, “Ya” atau “Tidak”).
    • Cabang (Branches): Jalur dari satu node ke node lain berdasarkan aturan keputusan.
  • Proses Pembelajaran: Algoritma Decision Tree memilih fitur dan nilai ambang (threshold) untuk membagi data sehingga memaksimalkan pemisahan kelas (untuk klasifikasi) atau meminimalkan kesalahan prediksi (untuk regresi). Kriteria pemilihan pembagian meliputi:
    • Gini Impurity: Mengukur ketidakmurnian kelas dalam node. Nilai Gini rendah menunjukkan node yang lebih “murni” (kebanyakan observasi dari satu kelas). \[ \text{Gini} = 1 - \sum_{i=1}^k p_i^2 \] di mana \(p_i\) adalah proporsi kelas \(i\).
    • Entropi/Information Gain: Mengukur pengurangan ketidakpastian setelah pembagian. Entropi dihitung sebagai: \[ \text{Entropy} = - \sum_{i=1}^k p_i \log_2(p_i) \] Information Gain adalah pengurangan entropi setelah pembagian.
    • Chi-Square (untuk data kategorikal): Dalam beberapa implementasi, uji chi-square digunakan untuk mengevaluasi signifikansi pembagian.
  • Contoh dalam Konteks Data Kategorikal: Misalkan kita memiliki data hubungan antara jenis kelamin dan kepercayaan terhadap kehidupan setelah kematian. Decision Tree dapat dibangun dengan:
    • Fitur: Jenis kelamin (pria/wanita).
    • Target: Kepercayaan (ya/tidak). Pohon mungkin membagi data berdasarkan jenis kelamin, lalu menentukan probabilitas kepercayaan untuk setiap kelompok.

Keunggulan Decision Tree

  • Interpretasi Mudah: Struktur pohon mudah dipahami, mirip dengan alur keputusan manusia.
  • Fleksibel: Dapat menangani data kategorikal dan kontinu tanpa memerlukan asumsi distribusi.
  • Menangani Interaksi: Secara otomatis menangkap interaksi antar variabel tanpa perlu menentukan istilah interaksi secara eksplisit.
  • Relevansi dengan Data Kategorikal: Cocok untuk data dalam tabel kontingensi, karena dapat memodelkan hubungan antar kategori.

Keterbatasan Decision Tree

  • Overfitting: Pohon yang terlalu dalam cenderung menangkap noise dalam data, menyebabkan kinerja buruk pada data baru.
  • Ketidakstabilan: Perubahan kecil dalam data dapat menghasilkan pohon yang sangat berbeda.
  • Bias terhadap Fitur Kategorikal: Fitur dengan banyak kategori cenderung dipilih lebih sering, meskipun tidak selalu relevan.
  • Tidak Cocok untuk Hubungan Linier: Kurang efektif jika hubungan antar variabel bersifat linier.

Pengendalian Overfitting

  • Pruning (Pemangkasan): Menghapus cabang yang tidak signifikan untuk menyederhanakan pohon.

  • Maksimum Kedalaman: Membatasi kedalaman pohon.

  • Minimum Sampel per Node: Memastikan setiap node memiliki jumlah observasi minimum.

  • Validasi Silang (Cross-Validation): Mengevaluasi kinerja model pada data yang tidak digunakan untuk pelatihan.

Random Forest

Random Forest adalah algoritma pembelajaran mesin berbasis ensemble yang menggabungkan banyak Decision Tree untuk meningkatkan akurasi dan stabilitas prediksi. Random Forest sangat populer untuk tugas klasifikasi dan regresi, termasuk dalam analisis data kategorikal, karena kemampuannya menangani data kompleks dan mengurangi overfitting.

Konsep Dasar Random Forest

  • Ensemble Learning: Random Forest membangun banyak pohon keputusan dan menggabungkan prediksinya melalui:
    • Mayoritas Voting (untuk klasifikasi): Kelas yang paling banyak dipilih oleh pohon-pohon menjadi prediksi akhir.
    • Rata-rata (untuk regresi): Nilai rata-rata dari semua pohon digunakan sebagai prediksi.
  • Dua Elemen Kunci Randomisasi:
    • Bootstrap Sampling (Bagging): Setiap pohon dilatih pada subset data yang diambil secara acak dengan pengembalian (bootstrap sample).
    • Pemilihan Fitur Acak: Pada setiap node, hanya subset acak dari fitur yang dipertimbangkan untuk pembagian.
  • Contoh dalam Konteks Data Kategorikal: Menggunakan data kepiting tapal kuda, Random Forest dapat memprediksi keberadaan satelit berdasarkan lebar karapas, warna, dan berat. Setiap pohon memproses subset data dan fitur, lalu hasilnya digabungkan untuk prediksi yang lebih akurat.

Keunggulan Random Forest

  • Akurasi Tinggi: Menggabungkan banyak pohon mengurangi varians dan meningkatkan kinerja prediksi dibandingkan Decision Tree tunggal.
  • Tahan terhadap Overfitting: Randomisasi dan bagging membuat model lebih robust terhadap noise.
  • Menangani Data Kompleks: Efektif untuk dataset besar, fitur campuran (kategorikal dan kontinu), dan hubungan non-linier.
  • Pentingnya Fitur: Random Forest dapat menghitung pentingnya fitur berdasarkan seberapa besar fitur tersebut mengurangi impurity (misalnya, Gini atau entropi) di seluruh pohon.
  • Estimasi Kesalahan OOB: Data out-of-bag menyediakan estimasi kinerja model tanpa memerlukan dataset validasi terpisah.
  • Relevansi dengan Data Kategorikal: Cocok untuk data seperti dalam tabel kontingensi, karena dapat menangani variabel kategorikal tanpa perlu transformasi ekstensif.

Keterbatasan Random Forest

  • Kurang Interpretabel: Berbeda dengan Decision Tree tunggal, Random Forest sulit diinterpretasikan karena melibatkan banyak pohon.
  • Komputasi Berat: Membutuhkan lebih banyak sumber daya komputasi untuk pelatihan dan prediksi, terutama untuk dataset besar.
  • Bias pada Fitur Kategorikal: Seperti Decision Tree, Random Forest dapat bias terhadap fitur dengan banyak kategori.
  • Tidak Cocok untuk Data Linear: Untuk hubungan linier sederhana, metode seperti regresi logistik sering lebih efisien.

Parameter yang Dapat Disesuaikan

  • Jumlah Pohon (n_estimators): Lebih banyak pohon meningkatkan akurasi tetapi menambah waktu komputasi.
  • Maksimum Fitur (max_features): Jumlah fitur yang dipertimbangkan pada setiap pembagian (biasanya \(\sqrt{\text{jumlah fitur}}\) untuk klasifikasi).
  • Kedalaman Pohon (max_depth): Membatasi kedalaman untuk mengurangi overfitting.
  • Minimum Sampel per Node: Mengontrol ukuran node untuk mencegah pohon terlalu spesifik.

Distribusi Probabilitas dalam Data Kategoris

Distribusi Bernoulli

Distribusi Bernoulli merupakan distribusi probabilitas diskrit untuk percobaan dengan dua hasil, yaitu sukses (1) dan gagal (0), dalam satu kali percobaan. Jika X menyatakan hasil dari suatu percobaan Bernoulli, maka X bernilai 1 dengan probabilitas p , dan bernilai 0 dengan probabilitas 1-p.

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

Keterangan Notasi:

  • \(X\): Variabel acak biner (bernilai 0 atau 1)
  • \(p\): Probabilitas terjadinya “sukses” (\(X = 1\)), dengan \(0 < p < 1\)
  • \(1 - p\): Probabilitas terjadinya “gagal” (\(X = 0\))

Distribusi Binomial

Distribusi binomial menggambarkan jumlah keberhasilan dalam n percobaan independen, di mana setiap percobaan memiliki probabilitas keberhasilan tetap (p) dan probabilitas kegagalan (q = 1 - p). Distribusi ini sangat penting dalam statistik karena sering digunakan dalam situasi dunia nyata, seperti uji coba kualitas produk, pemilihan acak, atau analisis kejadian biner.

Fungsi probabilitas dari distribusi binomial :

\[ P(X = x) = \binom{n}{x} p^x (1 - p)^{n - x}, \quad x = 0, 1, 2, \ldots, n \]

Keterangan Notasi:

  • \(X\): Variabel acak yang menyatakan jumlah keberhasilan
  • \(n\): Jumlah total percobaan
  • \(x\): Banyaknya keberhasilan yang diharapkan
  • \(p\): Probabilitas sukses dalam satu percobaan
  • \(\binom{n}{x}\): Kombinasi, dibaca “n choose x”, yaitu \(\displaystyle \binom{n}{x} = \frac{n!}{x!(n - x)!}\)

Syarat Distribusi Binomial

Agar suatu percobaan dapat dikategorikan sebagai distribusi binomial, harus memenuhi empat syarat berikut:

  1. Jumlah percobaan tetap (n): Jumlah percobaan (n) harus ditentukan sebelumnya.

  2. Dua hasil mungkin: Setiap percobaan hanya memiliki dua kemungkinan hasil, yaitu “sukses” atau “gagal”.

  3. Probabilitas konstan: Probabilitas keberhasilan (p) dalam setiap percobaan adalah konstan, sehingga probabilitas kegagalan adalah q = 1 - p.

  4. Independensi: Hasil dari satu percobaan tidak memengaruhi hasil percobaan lainnya.

Distribusi Multinomial

Distribusi multinomial menggambarkan probabilitas mendapatkan kombinasi tertentu dari hasil dalam n percobaan independen, di mana setiap percobaan dapat menghasilkan salah satu dari k kategori dengan probabilitas tetap untuk masing-masing kategori. Distribusi ini sering digunakan dalam situasi di mana ada lebih dari dua hasil yang mungkin, seperti dalam pemungutan suara, klasifikasi produk, atau analisis data kategorikal.

Fungsi probabilitas distribusi multinomial :

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

Dengan syarat:

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

Keterangan Notasi:

  • \(X_i\): Banyaknya kejadian kategori ke-\(i\)
  • \(x_i\): Realisasi dari \(X_i\), yaitu jumlah percobaan yang menghasilkan kategori ke-\(i\)
  • \(p_i\): Probabilitas terjadinya kategori ke-\(i\) dalam satu percobaan
  • \(n\): Jumlah total percobaan
  • \(k\): Jumlah kategori

Distribusi Poisson

Distribusi Poisson menggambarkan probabilitas terjadinya k kejadian dalam suatu interval waktu, ruang, atau unit pengamatan tertentu, di mana kejadian terjadi secara acak dan independen dengan rata-rata kejadian per interval (\(\lambda\)) yang konstan. Distribusi ini sering digunakan untuk memodelkan kejadian langka atau acak, seperti jumlah panggilan telepon dalam satu jam atau jumlah kecelakaan di suatu ruas jalan dalam sehari.

Fungsi probabilitas distribusi poisson :

\[ P(X = x) = (1 - p)^{x - 1} p, \quad x = 1, 2, 3, \ldots \]

Keterangan Notasi:

  • \(X\): Banyaknya percobaan sampai sukses pertama terjadi
  • \(x\): Realisasi dari \(X\) (keberhasilan pertama terjadi pada percobaan ke-\(x\))
  • \(p\): Probabilitas sukses dalam satu percobaan
  • \(1 - p\): Probabilitas gagal dalam satu percobaan

Syarat Distribusi Poisson

Agar suatu proses dapat dimodelkan dengan distribusi Poisson, harus memenuhi syarat berikut:

  1. Kejadian independen: Kejadian yang terjadi di satu interval tidak memengaruhi kejadian di interval lain.

  2. Laju rata-rata konstan: Rata-rata jumlah kejadian per unit waktu atau ruang (\(\lambda\)) tetap konstan.

  3. Kejadian langka: Probabilitas terjadinya lebih dari satu kejadian dalam interval waktu atau ruang yang sangat kecil mendekati nol.

  4. Interval terdefinisi: Pengamatan dilakukan dalam interval waktu, ruang, atau unit pengamatan yang jelas (misalnya, per jam, per kilometer, dll.).

Tabel Kontingensi 2x2

Tabel kontingensi 2x2 adalah salah satu alat dasar dalam statistik untuk menganalisis hubungan antara dua variabel kategorikal, masing-masing dengan dua kategori. Ini sering digunakan dalam penelitian survei, kesehatan, dan ilmu sosial untuk memahami apakah ada asosiasi antara variabel, seperti hubungan antara kebiasaan merokok dan penyakit tertentu.

Definisi dan Struktur

Tabel kontingensi 2x2 adalah matriks 2x2 yang menampilkan distribusi frekuensi dari dua variabel kategorikal biner. Struktur tabel ini terdiri dari:

Tabel kontingensi 2x2 memiliki struktur sebagai berikut :

\[ \begin{array}{|c|c|c|c|} \hline & \textbf{Kategori Positif} & \textbf{Kategori Negatif} & \textbf{Total} \\ \hline \textbf{Grup 1} & n_{11} & n_{12} & n_{1.} \\ \hline \textbf{Grup 2} & n_{21} & n_{22} & n_{2.} \\ \hline \textbf{Total} & n_{.1} & n_{.2} & n \\ \hline \end{array} \]

Keterangan Notasi:

Distribusi Peluang dalam Tabel Kontingensi 2x2

Peluang Bersama

Peluang bersama adalah probabilitas bahwa kedua variabel terjadi secara bersamaan dalam suatu sel tabel kontingensi:

\[ P(A_i, B_j) = \frac{n_{ij}}{n} \]

Peluang Marginal

Peluang marginal adalah probabilitas kejadian suatu variabel tanpa mempertimbangkan variabel lainnya : Peluang marginal baris :

\[ P(A_i) = \frac{n_{i.}}{n} \]

Peluang marginal kolom :

\[ P(B_j) = \frac{n_{.j}}{n} \]

Peluang Bersyarat

Peluang bersyarat adalah probabilitas suatu kejadian terjadi dengan syarat kejadian lain telah terjadi : \[ P(B_j|A_i) = \frac{P(A_i,B_j)}{P(A_i)} = \frac{n_{ij}}{n_i} \]

Jenis Kelamin Laki-Laki Perempuan Total
Laki-Laki 40 10 50
Perempuan 30 20 50
Total 70 30 100

Peluang Bersama

\[ P(\text{Laki-laki, Menyukai}) = \frac{40}{100} = 0.40 \]

\[ P(\text{Laki-laki, Tidak Menyukai}) = \frac{10}{100} = 0.10 \]

\[ P(\text{Perempuan, Menyukai}) = \frac{30}{100} = 0.30 \]

\[ P(\text{Perempuan, Tidak Menyukai}) = \frac{20}{100} = 0.20 \]

Peluang Marginal

\[ P(\text{Laki-laki}) = \frac{50}{100} = 0.50 \]

\[ P(\text{Perempuan}) = \frac{50}{100} = 0.50 \]

\[ P(\text{Menyukai}) = \frac{70}{100} = 0.70 \]

\[ P(\text{Tidak Menyukai}) = \frac{30}{100} = 0.30 \]

Peluang Bersyarat

\[ P(\text{Menyukai | Laki-laki}) = \frac{P(\text{Laki-laki, Menyukai})}{P(\text{Laki-laki})} = \frac{0.40}{0.50} = 0.80 \]

\[ P(\text{Menyukai | Perempuan}) = \frac{P(\text{Perempuan, Menyukai})}{P(\text{Perempuan})} = \frac{0.30}{0.50} = 0.60 \]

\[ P(\text{Tidak Menyukai | Laki-laki}) = \frac{P(\text{Laki-laki, Tidak Menyukai})}{P(\text{Laki-laki})} = \frac{0.10}{0.50} = 0.20 \]

\[ P(\text{Tidak Menyukai | Perempuan}) = \frac{P(\text{Perempuan, Tidak Menyukai})}{P(\text{Perempuan})} = \frac{0.20}{0.50} = 0.40 \]

Perhitungan dengan R :

# Data Observasi
library(lmerTest)
## Loading required package: lme4
## Loading required package: Matrix
## 
## Attaching package: 'lmerTest'
## The following object is masked from 'package:lme4':
## 
##     lmer
## The following object is masked from 'package:stats':
## 
##     step
data <- matrix(c(40,10,30,20), nrow = 2, byrow = TRUE)
colnames(data) <- c("Laki-Laki", "Perempuan")
rownames(data) <- c("Menyukai", "Tidak Menyukai")
n <- sum(data)

# Peluang Bersama
P_joint <- data / n

# Peluang Marginal
P_marginal_rows <- rowSums(data) / n
P_marginal_cols <- colSums(data) / n

# Peluang Bersyarat
P_conditional <- data / rowSums(data)

# Hasil
list(Peluang_Bersama = P_joint, Peluang_Marginal_Baris = P_marginal_rows, Peluang_Marginal_Kolom = P_marginal_cols,Peluang_Bersyarat = P_conditional)
## $Peluang_Bersama
##                Laki-Laki Perempuan
## Menyukai             0.4       0.1
## Tidak Menyukai       0.3       0.2
## 
## $Peluang_Marginal_Baris
##       Menyukai Tidak Menyukai 
##            0.5            0.5 
## 
## $Peluang_Marginal_Kolom
## Laki-Laki Perempuan 
##       0.7       0.3 
## 
## $Peluang_Bersyarat
##                Laki-Laki Perempuan
## Menyukai             0.8       0.2
## Tidak Menyukai       0.6       0.4

Interpretasi Hasil

  • Laki-laki cenderung lebih sering menyukai smartphone baru (40%) dibandingkan perempuan (30%), sementara proporsi perempuan yang tidak menyukai (20%) lebih tinggi dibandingkan laki-laki (10%).

  • Sampel terdiri dari jumlah laki-laki dan perempuan yang seimbang (masing-masing 50%), dan mayoritas responden (70%) menyukai smartphone baru, menunjukkan penerimaan produk yang cukup tinggi di antara kedua kelompok.

  • Laki-laki memiliki kecenderungan lebih tinggi untuk menyukai smartphone baru (80%) dibandingkan perempuan (60%). Ini menunjukkan bahwa jenis kelamin mungkin memengaruhi preferensi, dengan laki-laki lebih cenderung menyukai produk tersebut dibandingkan perempuan.

Ukuran Asosiasi dalam Data Kategori 2x2

Dalam analisis statistik, asosiasi dalam tabel kontingensi 2 × 2 digunakan untuk menentukan hubungan antara dua variabel kategori. Sebuah tabel kontingensi 2×2 menunjukkan frekuensi kejadian dua variabel dengan dua kategori, sehingga memungkinkan kita untuk menilai apakah terdapat hubungan statistik diantara mereka.

Risk Difference

Ukuran yang menunjukkan perbedaan mutlak dalam proporsi kejadian (misalnya, penyakit) antara dua kelompok, yaitu kelompok yang terpapar (misalnya, perokok) dan kelompok yang tidak terpapar (misalnya, bukan perokok).

Rumus dari Risk Difference :

\[ \text{RD} = \left( \frac{n_{11}}{n_{1.}} \right) - \left( \frac{n_{21}}{n_{2.}} \right) \]

Perhitungan dalam R :

RD <- function(n11, n12, n21, n22) {
(n11 / (n11 + n12)) - (n21 / (n21 + n22))
}
RD(40,10,30,20)
## [1] 0.2

Interpretasi : Laki-Laki menyukai smartphone 20% lebih tinggi dibandingkan perempuan menyukai smartphone

Relative RIsk

Ukuran yang membandingkan probabilitas suatu kejadian (misalnya, penyakit) pada kelompok yang terpapar suatu faktor risiko dengan probabilitas kejadian pada kelompok yang tidak terpapar. RR sering digunakan dalam studi kohort atau uji klinis untuk mengevaluasi apakah suatu paparan (seperti merokok) meningkatkan atau menurunkan risiko kejadian tertentu (seperti kanker paru-paru).

Rumus dari Relative Risk :

\[ \text{RD} = \frac{n_{11}}{n_{1.}} / \frac{n_{21}}{n_{2.}} \]

Perhitungan menggunakan R :

RR <- function(n11, n12, n21, n22) {
(n11 / (n11 + n12)) / (n21 / (n21 + n22))
}
RR(40,10,30,20)
## [1] 1.333333

Interpretasi : Laki-laki 1,33 kali lebih mungkin menyukai smartphone baru dibandingkan perempuan, menunjukkan adanya hubungan antara jenis kelamin dan preferensi.

Odds Ratio

Ukuran statistik yang membandingkan peluang (odds) terjadinya suatu kejadian (misalnya, penyakit) pada kelompok yang terpapar suatu faktor risiko dengan peluang kejadian pada kelompok yang tidak terpapar. OR sering digunakan dalam studi epidemiologi, khususnya studi kasus-kontrol, karena dapat diestimasi dari data retrospektif di mana risiko absolut tidak langsung tersedia.

Rumus Odds Ratio :

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

Perhitungan menggunakan R :

OR <- function(n11, n12, n21, n22) {
(n11 * n22) / (n12 * n21)
}
OR(40,10,30,20)
## [1] 2.666667

Interpretasi : Peluang laki-laki menyukai smartphone adalah sekitar 2,67 kali lebih besar dibandingkan perempuan.

Inferensi Tabel Kontingensi 2 Arah

Inferensi dalam statistik mengacu pada proses pengambilan kesimpulan mengenai populasi berdasarkan sampel data. Dalam konteks tabel kontingensi dua arah, inferensi digunakan untuk menganalisis hubungan antara dua variabel kategorikal yang disusun dalam tabel kontingensi.

Estimasi

Estimasi bertujuan untuk memperkirakan parameter populasi berdasarkan data sampel. Estimasi dibagi menjadi:

  • Estimasi Titik : Estimasi titik digunakan untuk menentukan satu nilai spesifk sebagai perkiraan terbaik dari parameter populasi.

  • Estimasi Interval : Estimasi interval bertujuan untuk memberikan rentang nilai yang diyakini mengandung parameter populasidengan tingkat kepercayaan tertentu.

Uji Hipotesis

Uji hipotesis adalah prosedur statistik untuk menguji pernyataan atau klaim tentang populasi berdasarkan data sampel. Dalam analisis data kategori, uji hipotesis sering digunakan untuk mengevaluasi hubungan antar variabel kategorik, seperti apakah ada asosiasi antara dua variabel atau apakah distribusi data sesuai dengan hasil yang diharapkan

Langkah Uji Hipotesis :

  1. Merumuskan Hipotesis Nol (H0) dan Hipotesis Alternatif (H1) H0: Pernyataan bahwa tidak ada efek atau hubungan (status quo). H1: Pernyataan bahwa ada efek atau hubungan.
  2. Menentukan Tingkat Signifikansi ($\alpha$) Biasanya $\alpha$ = 0.05, artinya kita menerima risiko 5% untuk menolak H0 yang sebenarnya benar.
  3. Menghitung Statistik Uji Statistik uji dihitung berdasarkan data dan metode yang dipilih (misalnya, chi-square).
  4. Menentukan p-value atau Nilai Kritis p-value menunjukkan probabilitas mendapatkan hasil setidaknya se-ekstrem yang diamati, jika H0 benar.
  5. Membuat Keputusan Jika p-value < $\alpha$, tolak H0. Jika p-value ≥ $\alpha$, gagal menolak H0.

Uji Proporsi

Uji proporsi adalah metode statistik untuk menguji hipotesis tentang proporsi populasi berdasarkan data sampel. Dalam analisis data kategori, uji ini sering digunakan untuk mengevaluasi apakah proporsi suatu kategori dalam populasi sesuai dengan nilai yang dihipotesiskan atau untuk membandingkan proporsi antar dua atau lebih kelompok.

Langkah Uji Proporsi :

  1. Merumuskan Hipotesis Nol (H0) dan Hipotesis Alternatif (H1) H0 : Proporsi populasi sama dengan nilai tertentu (misalnya, p = 0.5). H1: Proporsi populasi tidak sama dengan nilai tersebut (misalnya, p ≠ 0.5).
  2. Menentukan Tingkat Signifikansi ($\alpha$) Biasanya $\alpha$ = 0.05.
  3. Menghitung Statistik Uji Statistik : uji untuk proporsi biasanya adalah uji z.
  4. Menentukan p-value atau Nilai Kritis : p-value menunjukkan probabilitas mendapatkan hasil setidaknya se-ekstrem yang diamati, jika H0 benar.
  5. Membuat Keputusan : Jika p-value < $\alpha$, tolak H0. Jika p-value ≥ $\alpha$, gagal menolak H0.

Rumus statistik uji :

Untuk uji proporsi satu sampel, statistik uji z dihitung dengan rumus :

\[ z = \frac{p - p_0}{\sqrt{p_0(1 - p_0) \frac{1}{n}}} \]

Dimana :

  • \(p\) adalah proporsi sampel

  • \(p_0\) adalah proporsi yang dhipotesiskan dalam hipotesis nol

  • \(n\) adalah ukuran sampel

Asumsi :

Sampel diambil secara acak dan independen. Ukuran sampel cukup besar, biasanya \(np >= 5\) dan \(n(1-p) >= 5\)

Uji proporsi digunakan untuk membandingkan proporsi kejadian antara dua kelompok dalam tabel kontingensi, terutama untuk menentukan apakah terdapat perbedaan yang signifkan dalam proporsi kejadian antara dua kelompok yang berbeda

Contoh Soal :

Kondisi Sehat Tidak Sehat Total
Olahraga 30 10 40
Tidak Olahraga 20 40 60
Total 50 50 100

Hipotesis :

  • H0 : Tidak ada perbedaan proporsi antara 2 kelompok

  • H1 : Terdapat perbedaan proporsi antara 2 kelompok

Statistik Uji :

  1. Hitung Proporsi Sampel

\[ \hat{p_1} = \frac{30}{40} = 0.75 \] \[ \hat{p_2} = \frac{20}{60} = 0.333 \]

  1. Hitung proporsi gabungan

\[ \hat{p} = \frac{30 + 20}{40 + 60} = 0.50 \]

  1. Hitung statistik uji Z

\[ Z = \frac{0.75 - 0.333}{\sqrt{0.50(1 - 0.50)\left(\frac{1}{40} + \frac{1}{60}\right)}} \] \[ Z = \frac{0.417}{\sqrt{0.25(0.04167)}} = \frac{0.417}{0.2042} = 4.08 \]

Kesimpulan

Karena Z = 4.08 lebih besar dari nilai kritis untuk alpha = 0,05 (1.96) maka dapat disimpulkan bahwa menolak H0, yang berarti ada perbedaan signifikan antara 2 proporsi

Menggunakan perhitungan syntax R :

data <- matrix(c(30,10,20,40),nrow = 2, byrow = TRUE)
dimnames(data) <- list ("Olahraga" = c("Ya","Tidak"),"Sehat" = c("Ya","Tidak"))
print (data)
##         Sehat
## Olahraga Ya Tidak
##    Ya    30    10
##    Tidak 20    40
prop_test <- prop.test(x = c(data[1,1], data[2,1]),
n = c(sum(data[1,]), sum(data[2,])))
print(prop_test)
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  c(data[1, 1], data[2, 1]) out of c(sum(data[1, ]), sum(data[2, ]))
## X-squared = 15.042, df = 1, p-value = 0.0001052
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  0.2162937 0.6170396
## sample estimates:
##    prop 1    prop 2 
## 0.7500000 0.3333333

Kesimpulan

Karena p-value < 0,05, maka terdapat perbedaan proporsi antara kelompok sehat dan tidak sehat

Uji Asosiasi

Uji asosiasi adalah metode statistik yang digunakan untuk menentukan apakah ada hubungan yang signifikan secara statistik antara dua variabel kategorikal. Dalam konteks ini, kita akan fokus pada tabel kontingensi 2x2, yang merupakan bentuk sederhana dari tabel kontingensi. Uji ini sering digunakan dalam bidang epidemiologi, penelitian klinis, dan ilmu sosial untuk memahami apakah terdapat ketergantungan antara dua variabel, seperti hubungan antara merokok dan kanker paru-paru atau efektivitas pengobatan dan pemulihan.

Tujuan dari uji asosiasi adalah untuk mengukur kekuatan dan signifikansi hubungan antara dua variabel kategorikal. Dengan menggunakan uji ini, kita dapat menentukan apakah dua variabel tersebut saling berhubungan atau independen satu sama lain. Misalnya, dalam studi medis, uji asosiasi dapat membantu menentukan apakah paparan faktor risiko tertentu (seperti merokok) berkaitan dengan hasil tertentu (seperti penyakit).

  1. Risk Diference (RD) – Mengukur selisih risiko absolut antara dua kelompok

\[ RD = \left( \frac{n_{11}}{n_1} \right) - \left( \frac{n_{21}}{n_2} \right) \]

Standar Error :

\[ SE(RD) = \sqrt{\frac{\hat{p}_1(1 - \hat{p}_1)}{n_1} + \frac{\hat{p}_2(1 - \hat{p}_2)}{n_2}} \]

Statistik Uji Z :

\[ Z_{RD} = \frac{RD}{SE(RD)} \]

  1. Relative Risk (RR) – Mengukur perbandingan risiko antara dua kelompok.

\[ RR = \frac{n_{11}}{n_1} \, \Big/ \, \frac{n_{21}}{n_2} \]

Standard Error :

\[ SE(\ln RR) = \sqrt{\frac{1}{n_{11}} + \frac{1}{n_1} + \frac{1}{n_{21}} + \frac{1}{n_2}} \]

Statistik Uji Z :
\[ Z_{RR} = \frac{\ln RR}{SE(\ln RR)} \]

  1. Odds Ratio (OR) – Mengukur perbandingan odds antara dua kelompok

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

Standard Error :

\[ SE(\ln OR) = \sqrt{\frac{1}{n_{11}} + \frac{1}{n_{12}} + \frac{1}{n_{21}} + \frac{1}{n_{22}}} \]

Statistik Uji Z :

\[ Z_{OR} = \frac{\ln OR}{SE(\ln OR)} \]

Pengimplementasian pada contoh kasus :

Perhitungan RIsk Difference

\[ RD = \left( \frac{n_{11}}{n_1} \right) - \left( \frac{n_{21}}{n_2} \right) \]

\[ RD = \left( \frac{30}{40} \right) - \left( \frac{20}{60} \right) \]

\[ RD = 0.75 - 0.3333 = 0.4167 \]

Standard Error :

\[ SE(RD) = \sqrt{\frac{\hat{p}_1(1 - \hat{p}_1)}{n_1} + \frac{\hat{p}_2(1 - \hat{p}_2)}{n_2}} \]

\[ SE(RD) = \sqrt{\frac{0.75(1 - 0.75)}{40} + \frac{0.33(1 - 0.33)}{60}} \]

\[ SE(RD) = 0.09164 \]

Statistik Uji Z :

\[ Z_{RD} = \frac{RD}{SE(RD)} \]

\[ Z_{RD} = \frac{0.4167}{0.09164} \]

\[ Z_{RD} = 4.55 \]

Perhitungan Risk Relative

\[ RR = \frac{n_{11}}{n_1} \, \Big/ \, \frac{n_{21}}{n_2} \]

\[ RR = \frac{30}{40} \, \Big/ \, \frac{20}{60} \]

\[ RR = 2.25 \]

Standard Error :

\[ SE(\ln RR) = \sqrt{\frac{1}{n_{11}} + \frac{1}{n_1} + \frac{1}{n_{21}} + \frac{1}{n_2}} \]

\[ SE(\ln RR) = \sqrt{\frac{1}{30} + \frac{1}{40} + \frac{1}{20} + \frac{1}{60}} \]

\[ SE(\ln RR) = 0.2042 \]

Statistik Uji Z :

\[ Z_{RR} = \frac{\ln RR}{SE(\ln RR)} \]

\[ Z_{RR} = \frac{\ln 2.25}{0.2042} \]

\[ Z_{RR} = 2.88 \]

Perhitungan menggunakan syntax R

n11 <- 30; n12 <- 10; n21 <- 20; n22 <- 40
n1. <- n11 + n12; n2. <- n21 + n22

# Risk Difference
p1<-(n11/n1.)
p2<-(n21/n2.)
rd <- p1 - p2
se_rd <- sqrt((p1 * (1 - p1) / n1.) + p2*((1 - p2) / n2.))
z_rd <- rd / se_rd

# Relative Risk
rr <- (n11/n1.) / (n21/n2.)
se_ln_rr <- sqrt((1/n11) - (1/n1.) + (1/n21) - (1/n2.))
z_rr <- log(rr) / se_ln_rr

# Odds Ratio
or <- (n11 * n22) / (n12 * n21)
se_ln_or <- sqrt((1/n11) + (1/n12) + (1/n21) + (1/n22))
z_or <- log(or) / se_ln_or

# Hasil
list(RD = rd, SE_RD = se_rd, Z_RD = z_rd, RR = rr, SE_Ln_RR = se_ln_rr, Z_RR = z_rr, OR = or, SE_Ln_OR = se_ln_or, Z_OR = z_or)
## $RD
## [1] 0.4166667
## 
## $SE_RD
## [1] 0.09160351
## 
## $Z_RD
## [1] 4.548588
## 
## $RR
## [1] 2.25
## 
## $SE_Ln_RR
## [1] 0.2041241
## 
## $Z_RR
## [1] 3.97273
## 
## $OR
## [1] 6
## 
## $SE_Ln_OR
## [1] 0.4564355
## 
## $Z_OR
## [1] 3.925548

Uji Indepedensi

Uji independensi adalah metode statistik yang digunakan untuk menentukan apakah ada hubungan signifikan antara dua variabel kategorikal. Uji ini sering diterapkan pada tabel kontingensi, yang menampilkan distribusi frekuensi dari dua variabel. Tujuannya adalah untuk menguji hipotesis nol ((H_0)) bahwa kedua variabel independen (tidak ada hubungan) terhadap hipotesis alternatif ((H_1)) bahwa ada hubungan antara keduanya. Metode utama yang digunakan meliputi Uji Chi-Square, Uji Likelihood Ratio (G²), dan Uji Fisher’s Exact Test untuk sampel kecil.

Uji Chi-Square

Uji Chi-Square adalah metode yang paling umum untuk menguji independensi antara dua variabel kategorikal. Uji ini membandingkan frekuensi observasi (\(O_ij\)) dengan frekuensi yang diharapkan (\(E_ij\)) di bawah asumsi bahwa kedua variabel independen.

Rumus Statistik uji Chi-Square dihitung dengan :

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

Di mana:

\[ O_{ij}: \text{Frekuensi observasi pada sel baris } i \text{ dan kolom } j. \]

\[ E_{ij}: \text{Frekuensi yang diharapkan, dihitung dengan:} \]

\[ E_{ij} = \frac{n_i \times n_j}{n} \]

Derajat Bebas :

\[ df = (I-1) (J-1) \]

Interpretasi :

  • Nilai \(\chi^2\) dibandingkan dengan nilai kritis dari distribusi Chi-Square untuk \(df\) yang sesuai pada tingkat signifikansi \(\alpha\) (misalnya, \(\alpha = 0.05\) memberikan nilai kritis sekitar 3.841 untuk \(df = 1\)).

  • Jika \(\chi^2\) lebih besar dari nilai kritis atau p-value \(< \alpha\), tolak \(H_0\), menunjukkan bahwa variabel tidak independen (ada hubungan signifikan).

  • Jika p-value \(> \alpha\), gagal tolak \(H_0\), menunjukkan bahwa variabel kemungkinan independen.

Asumsi :

  • Data diambil secara acak dari populasi.

  • Frekuensi yang diharapkan (\(E_{ij}\)) pada setiap sel sebaiknya minimal 5 untuk tabel 2x2. Jika ada sel dengan \(E_{ij} < 5\), gunakan Uji Fisher’s Exact Test.

  • Tidak ada sel dengan frekuensi observasi nol.

Rumus Chi-Square :

\[ \chi^2 = \sum \frac{(O - E)^2}{E} \]

dimana :

  • \(O\) adalah nilai observasi dalam tabel kontingensi

  • \(E\) adalah nilai yang diharapkan, dihitung sebagai: \[ E_{ij} = \frac{(R_i \times C_j)}{N} \] dimana:

\(R_i\) = Total baris ke-\(i\)

\(C_j\) = Total kolom ke-\(j\)

\(N\) = Total sampel

Perhitungan manual pada tabel kontingensi 2x2 :

  1. Hitung nilai harapan (E) :

\[ E_{11} = \frac{(40 \times 50)}{100} = 20 \] \[ E_{12} = \frac{(40 \times 50)}{100} = 20 \] \[ E_{21} = \frac{(60 \times 50)}{100} = 30 \] \[ E_{22} = \frac{(60 \times 50)}{100} = 30 \]

  1. Hitung nilai Chi-Square :

\[ \chi^2 = \frac{(30 - 20)^2}{20} + \frac{(10 - 20)^2}{20} + \frac{(20 - 30)^2}{30} + \frac{(40 - 30)^2}{30} \]

\[ \chi^2 = \frac{100}{20} + \frac{100}{20} + \frac{100}{30} + \frac{100}{30} \]

\[ \chi^2 = 5 + 5 + \frac{10}{3} + \frac{10}{3} = 5 + 5 + 3.33 + 3.33 = 16.66 \]

Perhitungan menggunakan syntax R :

data <- matrix(c(30,10,20,40),nrow = 2, byrow = TRUE)
dimnames(data) <- list ("Olahraga" = c("Ya","Tidak"),"Kejadian" = c("Ya","Tidak"))
print (data)
##         Kejadian
## Olahraga Ya Tidak
##    Ya    30    10
##    Tidak 20    40
#Uji Chi-Square
chisq_test<- chisq.test(data)
print(chisq_test)
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  data
## X-squared = 15.042, df = 1, p-value = 0.0001052

Partisi Chi-Square

Partisi Chi-Square digunakan untuk mengidentifkasi kategori mana dalam tabel kontingensi yang bertanggung jawab atas hubungan yang signifkan. Jika uji Chi-Square pada tabel kontingensi I×J signifkan, maka partisi Chi-Square memungkinkan kita untuk menguraikan efek hubungan dalam subkelompok yang lebih kecil.

Contoh Kasus :

Pekerjaan Kantor Remote Hybrid Total
Teknik Mesin 120 80 100 300
Manajer 150 70 80 300
Staf 100 150 50 300
Total 370 300 230 900

Langkah 1 : Lakukan perhitungan uji Chi-Square

  1. Hitung nilai ekspektasinya

\[ E_{ij} = \frac{(n_i \times n_j)}{n} \]

  1. Hitung nilai Chi-Square

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

\[ \chi^2 = \frac{(120 - 123.33)^2}{123.33} + \frac{(150 - 123.33)^2}{123.33} + \frac{(100 - 123.33)^2}{123.33} + \frac{(150 - 100)^2}{100} + \frac{(70 - 100)^2}{100} + \frac{(80 - 76.67)^2}{76.67} + \frac{(100 - 76.67)^2}{76.67} + \frac{(150 - 76.67)^2}{76.67} + \frac{(80 - 76.67)^2}{76.67} \]

\[ \chi^2 = 2.77 + 7.48 + 6.50 + 25 + 9.00 + 0.11 + 9.38 + 51.02 + 0.11 = 111.37 \]

  1. Cari nilai dari derajat bebas

\[ (I - 1)(J - 1) = (3 - 1)(3 - 1) = 4 \]

Analisis Residual dalam Tabel Kontingensi

Residual dalam tabel kontingensi digunakan untuk mengidentifkasi sel mana yang menyumbang paling banyak terhadap hubungan antara variabel kategori. Residual mengukur selisih antara frekuensi yang diamati dan frekuensi yang diharapkan berdasarkan model independensi. Jika suatu sel memiliki residual yang besar (positif atau negatif), berarti frekuensi observasi dalam sel tersebut sangat berbeda dari yang diharapkan. Sebaliknya, jika residualnya kecil, berarti nilai observasi mendekati nilai ekspektasi, sehingga sel tersebut tidak banyak menyumbang terhadap hubungan antara variabel.

Jika residual 0 :

Tidak ada perbedaan signifkan antara jumlah observasi yang terjadi di suatu sel dengan jumlah yang diprediksi oleh model independensi. Artinya, variabel baris dan kolom dalam tabel tidak memiliki hubungan kuat atau tidak menunjukkan pola keterkaitan yang jelas.

  • Jika residual positif besar :

Frekuensi observasi jauh lebih tinggi dari yang diharapkan → menunjukkan hubungan positif yang kuat antara kategori tersebut.

  • Jika residual negatif besar :

Frekuensi observasi jauh lebih rendah dari yang diharapkan → menunjukkan hubungan negatifvatau tidak adanya asosiasi.

Covid (+) Covid (-) Total
Pakai Masker 450 200 650
Tidak Pakai Masker 50 300 350
Total 500 500 1000

Perhitungan menggunakan syntax R :

# Data Observasi
observed <- matrix(c(450, 200, 50, 300), nrow = 2, byrow = TRUE)
colnames(observed) <- c("Covid (+)", "Covid (-)")
rownames(observed) <- c("Pakai Masker", "Tidak Pakai Masker")

# Hitung nilai ekspektasi
expected <- chisq.test(observed)$expected

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

# Standardized Residual
row_sum <- rowSums(observed)
col_sum <- colSums(observed)
total_sum <- sum(observed)
standardized_residual <- (observed - expected) / sqrt(expected * (1 - row_sum / total_sum) * (1 - col_sum))
## Warning in sqrt(expected * (1 - row_sum/total_sum) * (1 - col_sum)): NaNs
## produced
# Menampilkan hasil
list(
Observed = observed,
Expected = expected,
Pearson_Residual = pearson_residual,
Standardized_Residual = standardized_residual
)
## $Observed
##                    Covid (+) Covid (-)
## Pakai Masker             450       200
## Tidak Pakai Masker        50       300
## 
## $Expected
##                    Covid (+) Covid (-)
## Pakai Masker             325       325
## Tidak Pakai Masker       175       175
## 
## $Pearson_Residual
##                    Covid (+) Covid (-)
## Pakai Masker        6.933752 -6.933752
## Tidak Pakai Masker -9.449112  9.449112
## 
## $Standardized_Residual
##                    Covid (+) Covid (-)
## Pakai Masker             NaN       NaN
## Tidak Pakai Masker       NaN       NaN

Tabel Kontingensi 3 Arah

Tabel kontingensi 3 arah adalah alat statistik yang digunakan untuk menganalisis hubungan antara tiga variabel kategorik. Tabel ini memperluas konsep tabel kontingensi 2 arah dengan menambahkan dimensi ketiga, sehingga memungkinkan analisis interaksi antar variabel. Tabel kontingensi 3 arah dijelaskan sebagai metode untuk mengevaluasi hubungan antar variabel dengan menggunakan distribusi frekuensi dan uji statistik seperti uji chi-square.

Tabel kontingensi 3 arah memiliki tiga dimensi, misalnya variabel \(A\), \(B\), dan \(C\), dengan kategori masing-masing \(I\), \(J\), dan \(K\). Frekuensi dalam sel tabel dinyatakan sebagai \(n_{ijk}\), di mana -i kategori untuk variabel \(A (1,2,3,...,I)\) kategrori untuk variabel \(B(1,2,3,...,J)\) dan \(k\) kategori untuk variabel \(C(1,2,3,...,K)\)

Tabel Parsial dan Marjinal

Struktur tabel kontingensi 3 arah :

Z = 1 Y = 1 Y = 2 Y = J Y = .
X = 1 n11 n121 n1J n1+.
X = 2 n12 n122 n1J+ n1++.
X = I n1I n121I n1JI n1I+.
Z = k
Z = K
X = . n++ n+++ n++++ .

Contoh Kasus :

  • Variabel A : Jenis kelamin (Pria,Wanita)

  • Variabel B : Status pekerjaan (Bekerja, Tidak Bekerja)

  • Variabel C : Tingkat pendidikan (Rendah, Tinggi)

\[ \begin{array}{|c|c|c|c|} \hline \textbf{Jenis Kelamin} & \textbf{Status Pekerjaan} & \textbf{Pendidikan Rendah} & \textbf{Pendidikan Tinggi} \\ \hline \text{Pria} & \text{Bekerja} & 50 & 30 \\ \text{Pria} & \text{Tidak Bekerja} & 20 & 10 \\ \text{Wanita} & \text{Bekerja} & 40 & 25 \\ \text{Wanita} & \text{Tidak Bekerja} & 15 & 10 \\ \hline \end{array} \]

Langkah penyelesaian :

  1. Hitung Frekuensi Marginak

    • Total Frekuensi (n) :

      \[ n = 50 + 30 + 20 + 10 + 40 + 25 + 15 + 10 = 200 \]

    • Marginal untuk Jenis kelamin (A) :

      \[ \text{Pria}: n_{1.1} = 50 + 30 + 20 + 10 = 110 \] \[ \text{Wanita}: n_{1.2} = 40 + 25 + 15 + 10 = 90 \]

    • Marginal untuk Status pekerjaan (B) :

      \[ \text{Bekerja}: n_{2.1} = 50 + 30 + 40 + 25 = 145 \] \[ \text{Tidak Bekerja}: n_{2.2} = 20 + 10 + 15 + 10 = 55 \]

    • Marginal untuk Tingkat pendidikan (C) :

    \[ \text{Rendah}: n_{3.1} = 50 + 20 + 40 + 15 = 125 \] \[ \text{Tinggi}: n_{3.2} = 30 + 10 + 25 + 10 = 75 \]

  2. Hitung Frekuensi Harapan (\(E_ijk\))

    Untuk setiap sel, gunakan rumus :

    \[ E_{ijk} = \frac{n_{i\cdot} \cdot n_{\cdot j} \cdot n_{\cdot k}}{n^2} \]Untuk sel (\(i = 1\), \(j = 1\), \(k = 1\)) (Pria, Bekerja, Pendidikan Rendah):

    \[ E_{111} = \frac{110 \cdot 145 \cdot 125}{200^2} = \frac{1,993,750}{40,000} = 49.84 \]

    Pria, Bekerja, Pendidikan Tinggi (\(E_{112}\)): \[ E_{112} = \frac{110 \cdot 145 \cdot 75}{200^2} = 29.89 \]

    Pria, Tidak Bekerja, Pendidikan Rendah (\(E_{121}\)): \[ E_{121} = \frac{110 \cdot 55 \cdot 125}{200^2} = 18.91 \]

    Pria, Tidak Bekerja, Pendidikan Tinggi (\(E_{122}\)): \[ E_{122} = \frac{110 \cdot 55 \cdot 75}{200^2} = 11.34 \]

    Wanita, Bekerja, Pendidikan Rendah (\(E_{211}\)): \[ E_{211} = \frac{90 \cdot 145 \cdot 125}{200^2} = 40.76 \]

    Wanita, Bekerja, Pendidikan Tinggi (\(E_{212}\)): \[ E_{212} = \frac{90 \cdot 145 \cdot 75}{200^2} = 24.46 \]

    Wanita, Tidak Bekerja, Pendidikan Rendah (\(E_{221}\)): \[ E_{221} = \frac{90 \cdot 55 \cdot 125}{200^2} = 15.49 \]

    Wanita, Tidak Bekerja, Pendidikan Tinggi (\(E_{222}\)): \[ E_{222} = \frac{90 \cdot 55 \cdot 75}{200^2} = 9.29 \]

    Sehingga dikumpulkan setiap hasil frekuensi ke dalam tabel berikut ini :

    \[ \begin{array}{|c|c|c|c|} \hline \text{Jenis Kelamin} & \text{Status Pekerjaan} & \text{Pendidikan Rendah} & \text{Pendidikan Tinggi} \\ \hline \text{Pria} & \text{Bekerja} & 49.84 & 29.89 \\ \text{Pria} & \text{Tidak Bekerja} & 18.91 & 11.34 \\ \text{Wanita} & \text{Bekerja} & 40.76 & 24.46 \\ \text{Wanita} & \text{Tidak Bekerja} & 15.49 & 9.29 \\ \hline \end{array} \]

    Perhitungan menggunakan syntax R :

    # Define the totals for rows, columns, and overall total
    n_total <- 200
    
    # Marginal totals for each dimension
    n_gender <- c(110, 90)  # Pria, Wanita
    n_employment <- c(145, 55)  # Bekerja, Tidak Bekerja
    n_education <- c(125, 75)  # Pendidikan Rendah, Pendidikan Tinggi
    
    # Calculate the expected frequencies for each combination
    E_111 <- (n_gender[1] * n_employment[1] * n_education[1]) / n_total^2
    E_112 <- (n_gender[1] * n_employment[1] * n_education[2]) / n_total^2
    E_121 <- (n_gender[1] * n_employment[2] * n_education[1]) / n_total^2
    E_122 <- (n_gender[1] * n_employment[2] * n_education[2]) / n_total^2
    E_211 <- (n_gender[2] * n_employment[1] * n_education[1]) / n_total^2
    E_212 <- (n_gender[2] * n_employment[1] * n_education[2]) / n_total^2
    E_221 <- (n_gender[2] * n_employment[2] * n_education[1]) / n_total^2
    E_222 <- (n_gender[2] * n_employment[2] * n_education[2]) / n_total^2
    
    # Store the expected frequencies in a matrix for easier display
    expected_frequencies <- matrix(c(E_111, E_112, E_121, E_122, E_211, E_212, E_221, E_222),
                                   nrow = 2, ncol = 2, byrow = TRUE,
                                   dimnames = list("Jenis Kelamin" = c("Pria", "Wanita"),
                                                   "Status Pekerjaan" = c("Bekerja", "Tidak Bekerja")))
    ## Warning in matrix(c(E_111, E_112, E_121, E_122, E_211, E_212, E_221, E_222), :
    ## data length differs from size of matrix: [8 != 2 x 2]
    # Print the matrix of expected frequencies
    print(expected_frequencies)
    ##              Status Pekerjaan
    ## Jenis Kelamin  Bekerja Tidak Bekerja
    ##        Pria   49.84375      29.90625
    ##        Wanita 18.90625      11.34375
  3. Hitung Statistik Chi-Square

    Dengan mengaplikasikan rumus dibawah ini :

    \[ \chi^2 = \sum \frac{(n_{ijk} - E_{ijk})^2}{E_{ijk}} \]

    \(E_{111}\) (Pria, Bekerja, Pendidikan Rendah): \[ \chi^2_{111} = \frac{(n_{111} - E_{111})^2}{E_{111}} = \frac{(50 - 49.84)^2}{49.84} = \frac{0.16^2}{49.84} = \frac{0.0256}{49.84} = 0.000514 \]

    \(E_{112}\) (Pria, Bekerja, Pendidikan Tinggi): \[ \chi^2_{112} = \frac{(n_{112} - E_{112})^2}{E_{112}} = \frac{(30 - 29.89)^2}{29.89} = \frac{0.11^2}{29.89} = \frac{0.0121}{29.89} = 0.000405 \]

    \(E_{121}\) (Pria, Tidak Bekerja, Pendidikan Rendah): \[ \chi^2_{121} = \frac{(n_{121} - E_{121})^2}{E_{121}} = \frac{(20 - 18.91)^2}{18.91} = \frac{1.09^2}{18.91} = \frac{1.1881}{18.91} = 0.0628 \]

    \(E_{122}\) (Pria, Tidak Bekerja, Pendidikan Tinggi): \[ \chi^2_{122} = \frac{(n_{122} - E_{122})^2}{E_{122}} = \frac{(10 - 11.34)^2}{11.34} = \frac{-1.34^2}{11.34} = \frac{1.7956}{11.34} = 0.1584 \]

    \(E_{211}\) (Wanita, Bekerja, Pendidikan Rendah): \[ \chi^2_{211} = \frac{(n_{211} - E_{211})^2}{E_{211}} = \frac{(40 - 40.76)^2}{40.76} = \frac{-0.76^2}{40.76} = \frac{0.5776}{40.76} = 0.0142 \]

    \(E_{212}\) (Wanita, Bekerja, Pendidikan Tinggi): \[ \chi^2_{212} = \frac{(n_{212} - E_{212})^2}{E_{212}} = \frac{(25 - 24.46)^2}{24.46} = \frac{0.54^2}{24.46} = \frac{0.2916}{24.46} = 0.0119 \]

    \(E_{221}\) (Wanita, Tidak Bekerja, Pendidikan Rendah): \[ \chi^2_{221} = \frac{(n_{221} - E_{221})^2}{E_{221}} = \frac{(15 - 15.49)^2}{15.49} = \frac{-0.49^2}{15.49} = \frac{0.2401}{15.49} = 0.0155 \]

    \(E_{222}\) (Wanita, Tidak Bekerja, Pendidikan Tinggi): \[ \chi^2_{222} = \frac{(n_{222} - E_{222})^2}{E_{222}} = \frac{(5 - 9.29)^2}{9.29} = \frac{-4.29^2}{9.29} = \frac{18.4041}{9.29} = 1.978 \]

    Perhitungan menggunakan syntax R :

# Observed frequencies (replace with actual data)
observed <- matrix(c(50, 30, 20, 10, 40, 25, 15, 5),
                   nrow = 2, ncol = 2, byrow = TRUE,
                   dimnames = list("Jenis Kelamin" = c("Pria", "Wanita"),
                                   "Status Pekerjaan" = c("Bekerja", "Tidak Bekerja")))
## Warning in matrix(c(50, 30, 20, 10, 40, 25, 15, 5), nrow = 2, ncol = 2, : data
## length differs from size of matrix: [8 != 2 x 2]
# Expected frequencies (calculated earlier)
expected <- matrix(c(49.84, 29.89, 18.91, 11.34, 40.76, 24.46, 15.49, 9.29),
                   nrow = 2, ncol = 2, byrow = TRUE,
                   dimnames = list("Jenis Kelamin" = c("Pria", "Wanita"),
                                   "Status Pekerjaan" = c("Bekerja", "Tidak Bekerja")))
## Warning in matrix(c(49.84, 29.89, 18.91, 11.34, 40.76, 24.46, 15.49, 9.29), :
## data length differs from size of matrix: [8 != 2 x 2]
# Calculate Chi-Square statistic for each cell
chi_square_values <- (observed - expected)^2 / expected

# Print the Chi-Square values for each cell
print(chi_square_values)
##              Status Pekerjaan
## Jenis Kelamin      Bekerja Tidak Bekerja
##        Pria   0.0005136437  0.0004048177
##        Wanita 0.0628291909  0.1583421517

Interpretasi :

  • Nilai Chi-Square yang lebih kecil menunjukkan bahwa data yang diamati lebih sesuai dengan model distribusi yang diharapkan, yang berarti tidak ada banyak perbedaan antara kategori.

  • Nilai Chi-Square yang lebih besar menunjukkan bahwa perbedaan antara data yang diamati dan yang diharapkan cukup signifikan, dan mungkin ada hubungan yang lebih kuat atau pola yang lebih jelas antara kategori-kategori tersebut.

Inferensi Tabel Kontingensi 3 Arah

Tabel kontingensi tiga arah digunakan untuk menganalisis hubungan antara dua variabel kategorik dengan mempertimbangkan variabel kontrol. Contohnya adalah hubungan antara kebiasaan merokok (X) dan kanker paru-paru (Y ) dengan variabel kontrol usia (Z).

Tabel ini terdiri dari beberapa tabel parsial (2×2) untuk setiap tingkat Z, serta tabel marginal yang mengabaikan Z. Ukuran asosiasi yang digunakan adalah odds ratio

Contoh Kasus :

Latar Belakang

Setiap perusahaan yang beroperasi 24 jam umumnya menerapkan pola kerja bergilir untuk menjaga proses produksi dan layanan. Penerapan sistem 2 shift dan 3 shift sering kali dipilih sesuai kebutuhan operasioanl, beban kerja serta jumlah tenaga kerja yang ada. Namun, perbedaan pola kerja gilir tersebut dapat memengaruhi kondisi fisik dan mental karyawan.

Jenis kelamin dan umur menjadi dua variabel utama, dimana karyawan berusia muda (< 30 tahun) biasanya lebih mudah beradaptasi secara fisik tetapi rentan pada aspek lain seperti kelelahan akibat kurangnya pengalaman. Sebaliknya, karyawan berusia lebih tua (> 30 tahun) biasanya memiliki ganguan keseharan lebih cepat ketika memiliki jam kerja yang tidak teratur. Selain faktor umur, jenis kelamin juga menjadi variabel pengaruh dimana fisik,mental dan tanggung jawab tiap gender memiliki proporsi yang berbeda

Untuk melihat bagaimana variabel umur dan jenis kelamin berpengaruh terhadap pilihan pola kerja gilir (2/3 shift) maka dibuat tabel kontingensi 3 arah. Melalui proses analisis tabel kontingensi 3 arah ini, perusahaan dapat mengevaluasi kebijakan penempatan shift, meminimalkan potensi kelelahan dan memastikan pemerataan beban kerja yang lebih proporsional dan adil

Tujuan Analisis

Tujuan analisis ini dilakukan sebagai berikut :

  • Mengetahui nilai asosiasi, dengan melihat nilai OR atau RR dapat menilai apakah pola kerja gilir 2 shfit memiliki kemungkinan lebih besar mengalami kelelahan dibandingkan kelompok kerja gilir 3 shift.

  • Membandingkan tingkat risiko antara kelompok, memperlihatkan terjadinya kelelahan antara kelompok pola kerja gilir 2 shift dengan 3 shift.

  • Menilai dampak terhadap keseluruhan, menunjukkan selisih risiko antara kelompok 2 shift dengan 3 shift

Alasan Confounder

Dalam analisis hubungan antara umur dan pola kerja gilir, jenis kelamin dapat berperan sebagai confounder, yaitu variabel yang memengaruhi hubungan kedua variabel utama. Jenis kelamin dapat memengaruhi pola kerja yang dipilih atau diberikan kepada karyawan, karena faktor-faktor seperti kemampuan fisik, kebijakan perusahaan, dan norma sosial. Misalnya, laki-laki mungkin lebih sering ditempatkan dalam pola kerja 2 shift karena dianggap memiliki daya tahan fisik yang lebih tinggi, sementara perempuan cenderung lebih banyak bekerja dalam pola 3 shift karena adanya kebijakan yang mempertimbangkan keseimbangan beban kerja dan faktor biologis. Selain itu, jenis kelamin juga dapat memengaruhi tingkat kelelahan yang dialami oleh karyawan dalam pola kerja tertentu. Perempuan, misalnya, mungkin menghadapi beban ganda antara pekerjaan dan tanggung jawab domestik yang dapat meningkatkan resiko kelelahan mereka dibandingkan dengan laki-laki

Proses Analisis

Jenis Kelamin Umur 2 Shift 3 Shift
Laki-Laki < 30 10 10
Laki-Laki > 30 39 40
Perempuan < 30 8 11
Perempuan > 30 26 39

Penulisan dalam syntax R :

data <- array(c(10,39,10,40,8,26,11,39),
 dim = c(2,2,2),
 dimnames = list(
 Umur = c("< 30 Tahun", "> 30 Tahun"),
 Pola_Kerja = c("2 Shift", "3 Shift"),
 Jenis_Kelamin = c("Laki-Laki","Perempuan")))

Perhitungan Parsial dan Marjinal

Tabel Parsial

Tabel parsial adalah tabel yang mengelompokkan X dan Y berdasarkan setiap level Z

#Tabel Parsial Berdasarkan Gender
freq_parsial_laki <- data[, , "Laki-Laki"]
freq_parsial_perempuan <- data[, , "Perempuan"]
#Tampilkan Hasil
freq_parsial_laki
##             Pola_Kerja
## Umur         2 Shift 3 Shift
##   < 30 Tahun      10      10
##   > 30 Tahun      39      40
freq_parsial_perempuan
##             Pola_Kerja
## Umur         2 Shift 3 Shift
##   < 30 Tahun       8      11
##   > 30 Tahun      26      39

Interpretasi tabel parsial :

Kelompok Laki-Laki :

  • Pada pola kerja 2 shift, laki-laki berusia > 30 tahun jauh lebih banyak dibandingkan laki-laki < 30 tahun

  • Pada pola kerja 3 shift, laki-laki berusia > 30 tahun lebih dominan dibandingkan laki-laki < 30 tahun

Kelompok Perempuan

  • Pada pola kerja 2 shift, perempuan berusia > 30 tahun lebih banyak dibandingkan perempuan < 30 tahun

  • Pada pola kerja 3 shift, perempuan berusia > 30 tahun lebih besar dibandingkan perempuan < 30 tahun

Kesimpulan Keseluruhan

Mayoritas karyawan di setiap kategori jenis kelamin yang bekerja pada pola kerja 2 maupun 3 shift adalah mereka yang berusia di atas 30 tahun.

Tabel Marjinal

Tabel marjinal adalah tabel yang mengabaikan variabel Z, dengan menjumlahkan data dari semua level z

freq_marjinal_x <- apply(data, 1, sum)
freq_marjinal_z <- apply(data, 3, sum)

#Tampilkan Hasil
freq_marjinal_x
## < 30 Tahun > 30 Tahun 
##         39        144
freq_marjinal_z
## Laki-Laki Perempuan 
##        99        84

Kesimpulan Tabel Marjinal

  • Pada tabel marjinal umur, nilai umur < 30 tahun lebih sedikit dibandingkan dengan karyawan umur > 30 tahun

  • Pada tabel marjinal jenis kelamin, karyawan dengan jenis kelamin laki laki lebih besar dibandingkan dengan perempuan

Peluang Bersama

Peluang bersama adalah probabilitas 2 atau lebih peristiwa yang terjadi secara bersama Rumus dari peluang bersama sebagai berikut :

\[ P(Z,X,Y)=f(Z,X,Y)/N \]

data <- data.frame(
 Jenis_Kelamin = rep(c("Laki-Laki", "Perempuan"), each = 2),
 Umur = c("< 30 Tahun", "> 30 Tahun"),
 dua_shift = c(10,39,8,26),
 total = c(20,79,19,65)
)
data$P_ZXY <- data$dua_shift / sum(data$total)
data

Kesimpulan Peluang Bersama

Kelompok Umur < 30 Tahun : Relatif sedikit responden yang berumur dibawah 30 tahun (masing - masing sekitar 5-6%) menunjukkan proporsi keseluruhan kelompok <30 memang lebih kecil dibandingkan > 30

Kelompok Umur > 30 Tahun : Responden diindikasikan bahwa mayoritas memang berada di kelompok > 30 tahun baik laki-laki maupun perempuan karena nilai peluang lebih besar

Perbandingan 2 Shift dan 3 Shift : Untuk umur > 30, diketahui bahwa proporsi laki laki di shift 2 hampir sama dengan shift 3 Lalu, untuk perempuan > 30 lebih banyak di 3 shift (0,2131) dibandingkan di 2 shift (0,1421)

Peluang Marjinal

Perhitungan peluang marjinal menggunakan syntax R :

# Gunakan data_cmh karena berbentuk array numerik

# Create data frame
data <- data.frame(
 Jenis_Kelamin = rep(c("Laki-Laki", "Perempuan"), each = 2),
 Umur = c("< 30 Tahun", "> 30 Tahun"),
 dua_shift = c(10,39,8,26),
 total = c(20,79,19,65)
)

# Select only numeric columns
numeric_data <- data[, sapply(data, is.numeric)]

# Compute total sum correctly
total <- sum(numeric_data)

# Compute joint probability matrix
joint_prob <- numeric_data / total

# Print results
print(joint_prob)
##    dua_shift      total
## 1 0.03759398 0.07518797
## 2 0.14661654 0.29699248
## 3 0.03007519 0.07142857
## 4 0.09774436 0.24436090
# Hitung peluang marjinal dengan benar
marginal_x <- apply(joint_prob, 1, sum)
marginal_z <- apply(joint_prob, 2, sum)

# Tampilkan hasil
marginal_x
## [1] 0.1127820 0.4436090 0.1015038 0.3421053
marginal_z
## dua_shift     total 
## 0.3120301 0.6879699

Kesimpulan

  • Dilihat dari umur, mayoritas sampel (78,96%) berumur > 30 tahun sementara hanya 21,31% yang < 30 tahun

  • Dilihat dari jenis kelamin, karyawan laki-laki (54,10%) sedikit lebih banyak daripada perempuan (45,90%)

Perhitungan Ukuran Asosiasi

Ukuran asosiasi dalam tabel kontingensi digunakan untuk mengukur kekuatan hubungan antara dua variabel kategori. Tiga ukuran asosiasi yang umum digunakan adalah Risk Diference (RD), Relative Risk (RR), dan Odds Ratio (OR).

Perhitungan menggunakan syntax R :

data <- matrix(c(10,10,8,11), nrow = 2, byrow = TRUE)
rownames(data) <- c("< 30 Tahun", "> 30 Tahun")
colnames(data) <- c("2 shift", "3 shift")

data
##            2 shift 3 shift
## < 30 Tahun      10      10
## > 30 Tahun       8      11

Risk Difference

Risk Difference (RD) adalah selisih antara risiko (probabilitas kejadian/outcome) pada kelompok yang terpapar (exposed) dan kelompok yang tidak terpapar (unexposed).

#Hitung Risk Difference
p1 <- data[1,1]/sum(data[1,])
p2 <- data[2,1]/sum(data[2,])
RD <- p1-p2
RD
## [1] 0.07894737

Interpretasi :

Berdasarkan nilai Risk Difference (RD) yang diperoleh, yaitu 0.07894737 (sekitar 7.89%), dapat disimpulkan bahwa risiko terjadinya kelelahan pada kelompok laki-laki pada shift 2 lebih tinggi 7.89% dibandingkan risiko pada perempuan pada shift 2.

Relative Risk

Relative Risk (RR) adalah ukuran asosiasi yang membandingkan risiko (probabilitas) terjadinya outcome pada kelompok terpapar (exposed) dengan risiko outcome pada kelompok tidak terpapar (unexposed).

RR <- p1/p2
RR
## [1] 1.1875

Interpretasi :

Relative Risk (RR) sebesar 1.1875 berarti bahwa kelompok laki-laki memiliki risiko terjadinya kelelahan yang 18.75% lebih tinggi dibandingkan kelompok perempuan.

Odds Ratio

Odds Ratio (OR) adalah ukuran asosiasi yang membandingkan odds (peluang) terjadinya suatu outcome (kejadian) di antara dua kelompok, umumnya antara kelompok terpapar (exposed) dan tidak terpapar (unexposed).

odds1 <- data[1,1]/data[1,2]
odds2 <- data[2,1]/data[2,2]
OR <- odds1/odds2
OR
## [1] 1.375

Interpretasi :

Odds Ratio (OR) = 1.375 menunjukkan bahwa odds (peluang) terjadinya kelelahan di kelompok laki- laki 1.375 kali lipat lebih besar daripada odds terjadinya outcome di kelompok perempuan

Indepedensi Bersyarat dalam Tabel Kontingensi 3 Arah

Independensi bersyarat terjadi ketika dua variabel dalam tabel kontingensi tidak memiliki pengaruh satu sama lain setelah memperhitungkan variabel ketiga yang digunakan sebagai kondisi atau pengontrol. Dalam konteks ini, kita mengatakan bahwa dua variabel independen secara bersyarat satu sama lain, mengingat variabel ketiga.

Untuk menguji independensi bersyarat dalam tabel kontingensi tiga arah, kita dapat menggunakan uji chi-square yang lebih kompleks, yaitu uji chi-square untuk tabel kontingensi tiga arah, atau pendekatan model statistik lainnya seperti log-linear models yang digunakan untuk menggambarkan hubungan antara lebih dari dua variabel kategorikal.

Pengujian Uji CMH

Uji Cochran-Mantel-Haenszel (Uji CMH) adalah uji statistik yang digunakan untuk menguji hubungan antara dua variabel kategorikal sambil mengontrol (atau memperhitungkan) satu atau lebih variabel pengganggu atau confounding.

Tujuan utama dari uji CMH adalah untuk menguji apakah ada hubungan antara dua variabel (biasanya variabel independen dan dependen) dalam kelompok yang berbeda, dengan mengontrol variabel ketiga yang diharapkan memengaruhi hubungan antara kedua variabel tersebut.

Uji CMH bekerja dengan memecah data menjadi strata atau lapisan berdasarkan variabel pengontrol (confounder). Setelah itu, uji chi-square dilakukan pada setiap lapisan untuk menguji apakah hubungan antara dua variabel utama tetap signifikan setelah memisahkan data berdasarkan variabel pengontrol.

Pengujian pada contoh kasus :

Hipotesis :

Hipotesis nol (H0) : \(\theta_{xy(k)} = 1\) untuk setiap \(k = 1, 2, \dots, K\)

Hipotesis satu (H1) : \(\theta_{xy(k)} \neq 1\) untuk paling sedikit satu \(k\)

Statistik Uji :

Statistik Cochran-Mantel-Haenszel (CMH) dirumuskan sebagai :

\[ CMH = \frac{\sum_k (n_{11k} - \mu_{11k})^2}{\sum_k \text{var}(n_{11k})} \]Proses perhitungan menggunakan syntax R :

data_cmh <- array(c(10,39,10,40,8,26,11,39),
 dim = c(2,2,2),
 dimnames = list(
 Umur = c("< 30 Tahun", "> 30 Tahun"),
 Pola_Kerja = c("2 Shift", "3 Shift"),
 Jenis_Kelamin = c("Laki-Laki","Perempuan")))

data_cmh
## , , Jenis_Kelamin = Laki-Laki
## 
##             Pola_Kerja
## Umur         2 Shift 3 Shift
##   < 30 Tahun      10      10
##   > 30 Tahun      39      40
## 
## , , Jenis_Kelamin = Perempuan
## 
##             Pola_Kerja
## Umur         2 Shift 3 Shift
##   < 30 Tahun       8      11
##   > 30 Tahun      26      39
cmh_base <- mantelhaen.test(data_cmh, correct = FALSE)
cmh_base
## 
##  Mantel-Haenszel chi-squared test without continuity correction
## 
## data:  data_cmh
## Mantel-Haenszel X-squared = 0.022132, df = 1, p-value = 0.8817
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
##  0.5175725 2.1541400
## sample estimates:
## common odds ratio 
##          1.055899

Interpretasi :

  • Karena p-value (0.8817) > 0.05, kita gagal menolak H0 Artinya, secara keseluruhan tidak ada bukti yang cukup untuk menyatakan bahwa common OR berbeda dari 1.

  • Common OR = 1.056 ≈ 1 menunjukkan bahwa, setelah mempertimbangkan strata (lapisan) yang dipisahkan odds terjadinya outcome di kelompok laki-laki dengan kelompok perempuan tidak jauh berbeda. Sehingga, asosiasi yang diamati antara variabel paparan dan outcome tidak signifikan secara statistik bila digabungkan lintas strata

Generalized Linear Model (GLM)

Generalized Linear Model (GLM) adalah perluasan dari regresi linier yang memungkinkan pemodelan data dengan distribusi non-normal dan hubungan non-linier antara variabel respons dan prediktor. GLM sangat fleksibel dan dapat digunakan untuk data dari keluarga eksponensial seperti binomial, Poisson, atau gamma, dengan fungsi tautan yang menghubungkan prediktor linier dengan rata-rata respons. Berdasarkan dokumen referensi, GLM sering digunakan dalam analisis data kesehatan, seperti data biner (sakit/tidak sakit), data hitungan (jumlah kasus penyakit), atau data kontinu dengan variansi tidak konstan.

Komponen Utama Generalized Linear Model (GLM)

Komponen Acak (Random Component)

Variabel respons \(Y\) mengikuti distribusi dari keluarga eksponensial, seperti normal, binomial, Poisson, atau gamma. Distribusi ini menentukan sifat data, misalnya biner untuk binomial atau hitungan untuk Poisson. Distribusi ini menentukan sifat data, misalnya biner untuk binomial atau hitungan untuk poisson.

Komponen Sistematik (Systematic Component)

Kombinasi linier dari variabel prediktor \(X_1, X_2, \dots, X_p\) yang dinyatakan sebagai:

\[ \eta = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \dots + \beta_p X_p \]

Di mana \(\eta\) adalah prediktor linier, dan \(\beta_0,\beta_1,\beta_2,.....\beta_p\) adalah koefisien regresi

Exponential Family

Exponential Family adalah kelompok distribusi probabilitas yang memiliki bentuk umum tertentu. Banyak distribusi yang sering digunakan dalam statistik, seperti distribusi normal, binomial, Poisson, dan gamma, termasuk dalam keluarga ini.

Beberapa ciri-ciri dari exponential family :

  1. Bentuk Umum: Fungsi distribusi probabilitas dapat diubah menjadi bentuk eksponensial, yang memudahkan analisis statistik, terutama dalam model-model statistik seperti regresi logistik dan model linear umum (GLM).
  1. Natural Parameter (\(\theta\)): Distribusi dalam exponential family selalu memiliki satu parameter yang disebut natural parameter, yang memudahkan kita untuk melakukan inferensi dan estimasi parameter melalui metode seperti maximum likelihood estimation (MLE).

  2. Fungsi Biasa (\(b(\theta)\)): Fungsi ini bergantung pada parameter \(\theta\) dan mengatur distribusi agar sesuai dengan bentuk eksponensial.

  3. Fungsi Normalisasi (\(c(y, \phi)\)): Fungsi ini mengatur agar distribusi tetap valid dengan cara memastikan bahwa total probabilitas dari semua kemungkinan nilai \(y\) adalah 1.

Distribusi termasuk dalam exponential family jika dapat ditulis dalam bentuk:

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

Contoh distribusi yang termasuk dalam exponential family:

  • Distribusi Normal

  • Distribusi Binomial

  • Distribusi Poisson

  • Distribusi Gamma

Contoh Pembuktian : Distribusi Binomial

Fungsi probabilitas distribusi binomial adalah:

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

Kita tuliskan ulang dalam bentuk exponential family:

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

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

Maka distribusi binomial termasuk dalam exponential family.

Model Regresi Logistik

Model regresi logistik adalah metode statistik yang digunakan untuk memodelkan hubungan antara variabel dependen kategorikal (biasanya biner, seperti 0 atau 1) dan satu atau lebih variabel independen. Berbeda dengan regresi linier yang memprediksi nilai kontinu, regresi logistik memprediksi probabilitas suatu peristiwa terjadi, seperti keberhasilan atau kegagalan, ya atau tidak, atau kategori lain yang bersifat diskrit. Model ini menggunakan fungsi logistik (atau sigmoid) untuk mengubah output linier menjadi probabilitas yang berkisar antara 0 dan 1. Menurut Hosmer dan Lemeshow (2013), regresi logistik sangat populer dalam berbagai bidang seperti kedokteran, pemasaran, dan ilmu sosial karena kemampuannya menangani data kategorikal dan interpretasi hasilnya yang intuitif.

Keunggulan regresi logistik adalah kemampuannya untuk menangani variabel independen yang bersifat kontinu maupun kategorikal, serta tidak memerlukan asumsi normalitas pada variabel independen seperti regresi linier. Namun, model ini tetap memiliki asumsi penting, seperti independensi observasi dan tidak adanya multikolinearitas antar variabel independen. Selain itu, interpretasi koefisien dalam bentuk odds ratio memudahkan pemahaman dampak variabel independen terhadap probabilitas kejadian. Sebagai contoh, odds ratio sebesar 2 untuk suatu variabel berarti bahwa peningkatan satu unit pada variabel tersebut meningkatkan odds kejadian sebanyak dua kali lipat, dengan asumsi variabel lain konstan (Hosmer & Lemeshow, 2013).

Model regresi logistik juga dapat diperluas menjadi regresi logistik multinomial untuk variabel dependen dengan lebih dari dua kategori atau regresi logistik ordinal untuk kategori yang memiliki urutan. Namun, tantangan seperti overfitting atau pelanggaran asumsi linearitas logit perlu diperhatikan. Untuk memastikan model yang baik, evaluasi seperti uji goodness-of-fit (misalnya, Hosmer-Lemeshow test) dan metrik akurasi seperti AUC-ROC sering digunakan. Dengan demikian, regresi logistik menjadi alat yang sangat berguna untuk analisis prediktif dalam konteks kategorikal.

Dalam regresi logistik, hubungan antara variabel independen dan probabilitas kejadian diungkapkan melalui persamaan logit, yaitu log dari odds ratio:

\[ \text{logit}(p) = \ln \left( \frac{p}{1-p} \right) = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \cdots + \beta_n x_n \]

di mana \(p\) adalah probabilitas kejadian, \(\beta_0\) adalah konstanta, dan \(\beta_1, \beta_2, \dots, \beta_n\) adalah koefisien regresi untuk variabel independen \(x_1, x_2, \dots, x_n\). Fungsi sigmoid,

\[ p = \frac{1}{1 + e^{-(\beta_0 + \beta_1 x_1 + \beta_2 x_2 + \cdots + \beta_n x_n)}} \]

memetakan nilai logit ke rentang probabilitas. Menurut Kleinbaum dan Klein (2010), estimasi parameter dalam model ini biasanya dilakukan menggunakan metode maximum likelihood estimation (MLE), yang mencari nilai koefisien yang memaksimalkan kemungkinan data yang diamati.

Contoh Penerapan Regresi Logistik

Regresi Logistik adalah aplikasi GLM untuk data biner, misalnya, memprediksi apakah pasien sakit (\(Y = 1\)) atau tidak sakit (\(Y = 0\)) berdasarkan usia (\(X_1\)) dan tekanan darah (\(X_2\)). Modelnya menggunakan fungsi tautan logit:

\[ \log \left( \frac{\mu}{1 - \mu} \right) = \beta_0 + \beta_1 X_1 + \beta_2 X_2 \]

Di mana \(\mu = P(Y = 1)\). Probabilitasnya adalah:

\[ \mu = \frac{e^{\beta_0 + \beta_1 X_1 + \beta_2 X_2}}{1 + e^{\beta_0 + \beta_1 X_1 + \beta_2 X_2}} \]

Simulasi data :

set.seed(123)
n <- 100
usia <- runif(n, 20, 80)
tekanan_darah <- runif(n, 90, 180)
X <- cbind(1, usia, tekanan_darah)
beta <- c(-5, 0.05, 0.02)
logit <- X %*% beta
prob <- exp(logit) / (1 + exp(logit))
y <- rbinom(n, 1, prob)
model <- glm(y ~ usia + tekanan_darah, family = binomial(link = "logit"))
summary(model)
## 
## Call:
## glm(formula = y ~ usia + tekanan_darah, family = binomial(link = "logit"))
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   -6.34827    1.75390  -3.620 0.000295 ***
## usia           0.05351    0.01473   3.634 0.000279 ***
## tekanan_darah  0.02940    0.01035   2.839 0.004522 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 137.19  on 99  degrees of freedom
## Residual deviance: 115.39  on 97  degrees of freedom
## AIC: 121.39
## 
## Number of Fisher Scoring iterations: 4

Interpretasi Koefisien Koefisien \(\beta_1\) dan \(\beta_2\) menunjukkan perubahan log-odds. Misalnya, jika \(\beta_1 = 0.05\), setiap kenaikan satu tahun usia meningkatkan log-odds sebesar 0.05. Odds ratio dihitung sebagai:}

\[ \text{Odds Ratio} = e^{\beta_1} \]

Jika $\(beta_1 = 0.05\) :

\[ e^{0.05} \approx 1.051 \]

Artinya, odds sakit meningkat 5.1 % per tahun usia.

Visualisasi Kurva Stigmoid

# Generate data
set.seed(123)
n <- 100
usia <- runif(n, 20, 80)  # Generating random age data between 20 and 80
tekanan_darah <- runif(n, 90, 180)  # Generating random blood pressure data between 90 and 180

# Combine usia and tekanan_darah into a data frame (X matrix)
X <- cbind(1, usia, tekanan_darah)

# Define beta coefficients
beta <- c(-5, 0.05, 0.02)

# Calculate logit (linear combination of inputs)
logit <- X %*% beta

# Calculate probability using the logistic function (sigmoid)
prob <- exp(logit) / (1 + exp(logit))

# Generate x values for visualization (evenly spaced 'Usia' values)
x_vals <- seq(min(usia), max(usia), length.out = 100)

# Calculate the logit values for the sigmoid curve
logit_vals <- beta[1] + beta[2] * x_vals + beta[3] * mean(tekanan_darah)

# Apply the logistic function to get the probability values for the sigmoid curve
prob_vals <- exp(logit_vals) / (1 + exp(logit_vals))

# Plot the sigmoid curve
plot(x_vals, prob_vals, type = "l", col = "blue", lwd = 2, 
     xlab = "Usia", ylab = "Probability", main = "Sigmoid Curve Visualization")
grid()  # Add grid lines to the plot

Contoh 2 Penerapan Regresi Logistik

Jumlah kecelakaan lalu lintas di 50 kota selama satu bulan. Variabel-variabelnya adalah:

  • Kecelakaan (\(Y\)): Jumlah kecelakaan lalu lintas (variabel respon, count data).

  • Kepadatan \((X_1\)): Kepadatan penduduk kota (dalam ribu jiwa per km², variabel kontinu).

  • Curah Hujan (\(X_2\)): Rata-rata curah hujan bulanan (dalam mm, variabel kontinu).

# Set seed untuk reproduksibilitas
set.seed(123)
n <- 50
data_kecelakaan <- data.frame(
 Kecelakaan = rpois(n, lambda = 10 + 2 * runif(n, 1, 5)), # Simulasi count data
 Kepadatan = runif(n, 1, 10), # Kepadatan: 1-10 ribu/km²
 CurahHujan = runif(n, 50, 200) # Curah hujan: 50-200 mm
)
# Lihat beberapa baris data
head(data_kecelakaan)
# Ringkasan data
summary(data_kecelakaan)
##    Kecelakaan      Kepadatan       CurahHujan    
##  Min.   : 4.00   Min.   :1.429   Min.   : 57.93  
##  1st Qu.:13.00   1st Qu.:3.499   1st Qu.:104.41  
##  Median :15.00   Median :5.702   Median :133.03  
##  Mean   :15.36   Mean   :5.666   Mean   :133.27  
##  3rd Qu.:18.00   3rd Qu.:7.248   3rd Qu.:165.66  
##  Max.   :26.00   Max.   :9.871   Max.   :195.40
# Visualisasi distribusi variabel respon
hist(data_kecelakaan$Kecelakaan, main = "Distribusi Jumlah Kecelakaan",
 xlab = "Jumlah Kecelakaan", col = "lightblue", breaks = 15)

Eksplorasi data dengan Visualisasi :

# Scatter plot: Kecelakaan vs Kepadatan
plot(data_kecelakaan$Kepadatan, data_kecelakaan$Kecelakaan,
 main = "Kecelakaan vs Kepadatan Penduduk",
 xlab = "Kepadatan (ribu/km²)", ylab = "Jumlah Kecelakaan",
 pch = 16, col = "blue")

# Scatter plot: Kecelakaan vs Curah Hujan
plot(data_kecelakaan$CurahHujan, data_kecelakaan$Kecelakaan,
 main = "Kecelakaan vs Curah Hujan",
 xlab = "Curah Hujan (mm)", ylab = "Jumlah Kecelakaan",
 pch = 16, col = "green")

Scatter plot membantu mengidentifikasi apakah ada pola linier atau non-linier antara prediktor dan respon.

Melakukan Pemodelan GLM

model_glm <- glm(Kecelakaan ~ Kepadatan + CurahHujan,
 family = poisson(link = "log"),
 data = data_kecelakaan)
# Ringkasan model
summary(model_glm)
## 
## Call:
## glm(formula = Kecelakaan ~ Kepadatan + CurahHujan, family = poisson(link = "log"), 
##     data = data_kecelakaan)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) 2.4807660  0.1673138  14.827   <2e-16 ***
## Kepadatan   0.0155822  0.0160207   0.973    0.331    
## CurahHujan  0.0012092  0.0009221   1.311    0.190    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 60.644  on 49  degrees of freedom
## Residual deviance: 58.246  on 47  degrees of freedom
## AIC: 291.09
## 
## Number of Fisher Scoring iterations: 4

Model: \[ \log(\mu_i) = \beta_0 + \beta_1 \text{Kepadatan}_i + \beta_2 \text{CurahHujan}_i \]

Fungsi penghubung log: \[ \mu_i = \exp(\beta_0 + \beta_1 \text{Kepadatan}_i + \beta_2 \text{CurahHujan}_i) \]

Interpretasi Koefisien :

Jika \(\beta_1 > 0\), kenaikan 1 unit pada Kepadatan meningkatkan \(\log(\mu)\) sebesar \(\beta_1\), atau \(\mu\) meningkat secara multiplikatif sebesar \(\exp(\beta_1)\).

P-value \(<\) \(\alpha\) (0.05) menunjukkan prediktor signifikan.

Melakukan Evaluasi Model

# Deviasi model
deviance(model_glm)
## [1] 58.24589
df.residual(model_glm)
## [1] 47

Kesimpulan :

  • Deviasi Model (58.24589) menunjukkan bahwa model yang digunakan memiliki deviasi yang cukup rendah, yang mengindikasikan bahwa model tersebut mampu mencocokkan data yang diamati dengan baik. Namun, untuk memastikan kualitas model secara keseluruhan, perlu dilakukan perbandingan dengan model lain menggunakan metrik yang sama atau tambahan lainnya.

  • Derajat Kebebasan Residual (47) menunjukkan jumlah data yang digunakan untuk memperkirakan residual model. Nilai ini menunjukkan bahwa model memiliki derajat kebebasan yang wajar, yang berarti model memiliki cukup data untuk menghitung kesalahan residual secara akurat.

Melakukan Plot Residual

# Plot residual
plot(fitted(model_glm), residuals(model_glm, type = "deviance"),
 main = "Residual vs Fitted Values",
 xlab = "Nilai Prediksi", ylab = "Residual Deviance",
 pch = 16, col = "red")
abline(h = 0, lty = 2)

Melakukan Prediksi Nilai

# Prediksi nilai
data_kecelakaan$Prediksi <- predict(model_glm, type = "response")
# Plot efek Kepadatan
plot(data_kecelakaan$Kepadatan, data_kecelakaan$Kecelakaan,
 main = "Kecelakaan: Observasi vs Prediksi",
 xlab = "Kepadatan (ribu/km²)", ylab = "Jumlah Kecelakaan",
 pch = 16, col = "blue")
points(data_kecelakaan$Kepadatan, data_kecelakaan$Prediksi,
 col = "red", pch = 16)
legend("topleft", legend = c("Observasi", "Prediksi"),
 col = c("blue", "red"), pch = 16)

# Plot efek Curah Hujan
plot(data_kecelakaan$CurahHujan, data_kecelakaan$Kecelakaan,
 main = "Kecelakaan: Observasi vs Prediksi",
 xlab = "Curah Hujan (mm)", ylab = "Jumlah Kecelakaan",
 pch = 16, col = "green")
points(data_kecelakaan$CurahHujan, data_kecelakaan$Prediksi,
 col = "red", pch = 16)
legend("topleft", legend = c("Observasi", "Prediksi"),
 col = c("green", "red"), pch = 16)

Model Regresi Poisson

Regresi Poisson adalah metode statistik yang digunakan untuk menganalisis data hitung, yaitu data berupa bilangan bulat non-negatif yang menunjukkan frekuensi kejadian, seperti jumlah pelanggan yang datang per jam atau jumlah klik pada situs web per hari. Model ini mengasumsikan bahwa data mengikuti distribusi Poisson, di mana rata-rata dan variansnya sama. Hal ini membuat regresi Poisson ideal untuk data yang tidak negatif dan sering kali memiliki distribusi miring ke kanan, yang umum ditemukan dalam bidang seperti epidemiologi, ekonomi, dan ilmu sosial. Berbeda dengan regresi kuadrat terkecil biasa yang mengasumsikan distribusi normal dan varians konstan, regresi Poisson lebih sesuai untuk data hitung dengan varians yang meningkat seiring rata-rata.

Regresi Poisson memodelkan logaritma natural dari jumlah harapan (expected count) sebagai kombinasi linier dari variabel-variabel prediktor, memastikan bahwa prediksi tetap positif, yang penting untuk data hitung. Secara matematis, untuk observasi ke-\(i\), model dinyatakan sebagai:

\[ \log(\lambda_i) = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_p x_{ip} \]

di mana \(\lambda_i\) adalah jumlah harapan, \(\beta_0\) adalah intersep, dan \(\beta_1, \beta_2, \dots, \beta_p\) adalah koefisien untuk variabel prediktor \(x_{i1}, x_{i2}, \dots, x_{ip}\). Dengan mengekspansikan kedua sisi, jumlah harapan diperoleh:

\[ \lambda_i = e^{\beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_p x_{ip}} \]

Koefisien \(\beta_k\) menunjukkan perubahan pada log-jumlah harapan untuk setiap kenaikan satu unit pada \(x_{ik}\). Dengan mengekspansikan \(\beta_k\), kita mendapatkan incidence rate ratio, yang menunjukkan faktor perkalian jumlah harapan untuk kenaikan satu unit pada variabel prediktor. Misalnya, jika \(\beta_1 = 0.01798\) untuk jam sejak mesin dibersihkan, jumlah cacat meningkat sekitar 1.8% per jam (\(e^{0.01798} \approx 1.018\)), signifikan pada p-value kurang dari 0.05.

Asumsi dan Keterbatasan

  1. Respons Poisson : Hitungan per satuan waktu atau ruang dan mengikuti disribusi poisson
  2. Indepedensi : Observasi bebas, artinya satu obeservasi tidak mememngaruhi yang lain
  3. Rata-rata = Varians : Rata - rata dan varians harus sama pada setiap observasi
  4. Linearitas : Logaritma jumlah harapan merupakan fungsi linear dari variabel prediktor

Namun, pelanggaran asumsi sering terjadi, terutama overdispersi, di mana varians melebihi rata-rata, sehingga model Poisson tidak lagi sesuai. Dalam kasus ini, alternatif seperti regresi binomial negatif atau regresi kuasi-Poisson (dengan parameter dispersi, misalnya \((ϕ=4.447)\) dapat digunakan. Overdispersi dapat menyebabkan standar eror yang terlalu kecil, sehingga koefisien mungkin tidak signifikan pada taraf 0,05. Selain itu, jika data memiliki terlalu banyak nilai nol (excess zeros), model Poisson dengan inflasi nol (zero-inflated Poisson, ZIP) dapat digunakan, yang memisahkan pemodelan proses hitung dan probabilitas nol berlebih

Contoh Penerapan Regresi Poisson

Untuk data hitungan, seperti jumlah kasus penyakit di suatu wilayah, digunakan regresi Poisson dengan fungsi tautan log :

\[ \log(\mu) = \beta_0 + \beta_1 X_1 + \beta_2 X_2 \]

Di mana \(miu\) adalah rata-rata jumlah kasus. Untuk mengilustrasikan penggunaan offset, kita asumsikan \(X_1\) adalah populasi (sebagai exposure) dan \(X_2\) adalah tingkat polusi.

set.seed(456)
n <- 100
populasi <- runif(n, 1000, 5000)
polusi <- runif(n, 0, 10)
X <- cbind(1, polusi)
beta <- c(1, 0.1)
log_mu <- log(populasi) + X %*% beta
mu <- exp(log_mu)
y <- rpois(n, mu)

# Correctly combine y and polusi (or X)
data <- data.frame(y, polusi)  # or you can use X if you'd like to include the intercept too
model_poisson <- glm(y ~ polusi + offset(log(populasi)), family = poisson(link = "log"))
summary(model_poisson)
## 
## Call:
## glm(formula = y ~ polusi + offset(log(populasi)), family = poisson(link = "log"))
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept) 1.0018489  0.0020228   495.3   <2e-16 ***
## polusi      0.0996899  0.0002896   344.3   <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: 124240.731  on 99  degrees of freedom
## Residual deviance:     87.872  on 98  degrees of freedom
## AIC: 1228.4
## 
## Number of Fisher Scoring iterations: 3

Diagnostik dan Overdispersion

Untuk memvalidasi model, periksa: \[ Deviance: \quad D = 2 \sum_{i=1}^{n} \left[ \ell(y_i, y_i) - \ell(\mu_i, y_i) \right] \]

\[ Residual: \quad r_i = \frac{y_i - \hat{\mu}_i}{\sqrt{\hat{\mu}_i}} \]

\[ Overdispersion: \quad \phi = \frac{\sum (y_i - \hat{\mu}_i)^2}{n - p} \]

deviance_poisson <- deviance(model_poisson)
df_residual <- df.residual(model_poisson)
dispersion <- sum((resid(model_poisson, type = "pearson")^2)) / df_residual
list(deviance = deviance_poisson, dispersion = dispersion)
## $deviance
## [1] 87.8722
## 
## $dispersion
## [1] 0.896745

Interpretasi :

Karena hasil dispersion < 1 dimana 0.896745 < 1 maka overdispersion tidak terjadi

Contoh 2 Pengerjaan Regresi Poisson

# Set seed untuk reproduksibilitas
set.seed(456)
# Buat dataset fiktif
n <- 60
data_toko <- data.frame(
 Pengunjung = rpois(n, lambda = 20 + 5 * sample(0:3, n, replace = TRUE) + 0.05 * runif(n, 100, 500)),
 Promosi = sample(0:3, n, replace = TRUE),
 UkuranToko = runif(n, 100, 500)
)
# Eksplorasi data
head(data_toko)
summary(data_toko)
##    Pengunjung       Promosi        UkuranToko   
##  Min.   :25.00   Min.   :0.000   Min.   :104.1  
##  1st Qu.:36.75   1st Qu.:0.000   1st Qu.:197.7  
##  Median :45.00   Median :1.000   Median :308.2  
##  Mean   :44.65   Mean   :1.367   Mean   :303.1  
##  3rd Qu.:51.25   3rd Qu.:2.000   3rd Qu.:421.0  
##  Max.   :68.00   Max.   :3.000   Max.   :495.0

Berdasarkan hasil summary data yang ditampilkan, dapat disimpulkan bahwa dataset ini menggambarkan variabel yang berkaitan dengan jumlah pengunjung, promosi, dan ukuran toko yang memiliki karakteristik distribusi yang berbeda. Berikut adalah beberapa hal yang dapat disimpulkan secara keseluruhan:

  1. Pengunjung:
-   Distribusi data pengunjung relatif merata, dengan rata-rata sekitar 44.65 dan median di angka 45.00. Hal ini menunjukkan bahwa jumlah pengunjung cenderung seimbang di sekitar nilai tengah dan tidak terlalu terfokus pada angka tertentu.

-   Penyebaran data yang cukup luas terlihat dengan adanya variasi antara nilai minimum dan maksimum yang cukup jauh.
  1. Promosi:
-   Variabel promosi memiliki rentang nilai yang lebih sempit dengan dominasi nilai 0 dan 1. Hal ini dapat diartikan bahwa sebagian besar toko dalam dataset ini tidak memberikan promosi yang signifikan atau hanya memberikan sedikit promosi.

-   Rata-rata yang cukup rendah (1.367) dan median yang berada di angka 1 menunjukkan bahwa promosi tidak berperan besar dalam banyak data yang tercatat.
  1. UkuranToko:

    • Ukuran toko memiliki rentang yang lebih luas, dengan nilai rata-rata 303.1 dan median 308.2. Ini menunjukkan bahwa meskipun mayoritas toko memiliki ukuran yang serupa, terdapat beberapa toko dengan ukuran yang lebih besar.

    • Variasi ukuran toko yang cukup besar (dengan nilai minimum 104.1 dan maksimum 495.0) menunjukkan adanya perbedaan signifikan dalam kapasitas atau luas toko yang ada dalam dataset ini.

hist(data_toko$Pengunjung, main = "Distribusi Jumlah Pengunjung Harian",
 xlab = "Jumlah Pengunjung", col = "lightgreen", breaks = 15)

# Visualisasi eksplorasi
par(mfrow = c(1, 2))
boxplot(Pengunjung ~ Promosi, data = data_toko,
 main = "Pengunjung vs Jumlah Promosi",
 xlab = "Jumlah Promosi", ylab = "Jumlah Pengunjung",
 col = "lightblue")
plot(data_toko$UkuranToko, data_toko$Pengunjung,
 main = "Pengunjung vs Ukuran Toko",
 xlab = "Ukuran Toko (m²)", ylab = "Jumlah Pengunjung",
 pch = 16, col = "purple")

Interpretasi :

  1. Boxplot: Pengunjung vs Jumlah Promosi
    • Grafik ini menunjukkan distribusi jumlah pengunjung berdasarkan jumlah promosi yang dilakukan (0, 1, 2, dan 3).

    • Terdapat perbedaan yang cukup kecil pada distribusi jumlah pengunjung untuk berbagai jumlah promosi. Hasil ini mengindikasikan bahwa jumlah promosi tidak memiliki pengaruh yang signifikan terhadap jumlah pengunjung.

    • Sebagai catatan, nilai pengunjung yang paling rendah, tertinggi, dan rata-rata hampir seragam di seluruh kategori jumlah promosi (0, 1, 2, dan 3).

  2. Scatter Plot: Pengunjung vs Ukuran Toko
    • Grafik ini menunjukkan hubungan antara jumlah pengunjung dan ukuran toko (dalam satuan meter persegi).

    • Tampak bahwa hubungan antara kedua variabel ini tidak terlalu jelas atau linear. Ada distribusi yang tersebar dengan tidak ada pola yang jelas. Ini mengindikasikan bahwa ukuran toko mungkin tidak memiliki pengaruh langsung yang signifikan terhadap jumlah pengunjung.

par(mfrow = c(1, 1))
# Pemasangan model
model_poisson <- glm(Pengunjung ~ Promosi + UkuranToko,
 family = poisson(link = "log"),
 data = data_toko)
summary(model_poisson)
## 
## Call:
## glm(formula = Pengunjung ~ Promosi + UkuranToko, family = poisson(link = "log"), 
##     data = data_toko)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  3.8804531  0.0567285  68.404   <2e-16 ***
## Promosi      0.0263013  0.0168861   1.558   0.1193    
## UkuranToko  -0.0003925  0.0001646  -2.384   0.0171 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 141.15  on 59  degrees of freedom
## Residual deviance: 133.41  on 57  degrees of freedom
## AIC: 476.22
## 
## Number of Fisher Scoring iterations: 4

Interpretasi :

Model ini menunjukkan bahwa jumlah promosi tidak memiliki pengaruh signifikan terhadap jumlah pengunjung, dengan p-value sebesar 0.1193 yang lebih besar dari 0.05. Hal ini mengindikasikan bahwa peningkatan jumlah promosi tidak berhubungan langsung dengan perubahan jumlah pengunjung yang signifikan. Di sisi lain, ukuran toko terbukti memiliki pengaruh signifikan terhadap jumlah pengunjung, dengan p-value 0.0171 yang lebih kecil dari 0.05. Koefisien negatif pada ukuran toko (-0.0003925) menunjukkan bahwa setiap penambahan satu meter persegi ukuran toko berhubungan dengan penurunan jumlah pengunjung sekitar 0.039%, yang menunjukkan bahwa toko yang lebih besar cenderung menarik lebih sedikit pengunjung.

Secara keseluruhan, meskipun jumlah promosi tidak menunjukkan pengaruh yang signifikan, ukuran toko memiliki pengaruh yang signifikan terhadap jumlah pengunjung. Model ini juga menunjukkan penurunan deviance dari 141.15 (null deviance) menjadi 133.41 (residual deviance), yang menunjukkan bahwa model memberikan penurunan yang signifikan dalam menjelaskan variabilitas data. Nilai AIC yang sebesar 476.22 menunjukkan bahwa model ini cukup baik dalam memprediksi jumlah pengunjung berdasarkan variabel yang digunakan.

Evaluasi dan Overdispersion

# Evaluasi model
deviance(model_poisson)
## [1] 133.4076
dispersion <- deviance(model_poisson) / df.residual(model_poisson)
cat("Dispersion parameter:", dispersion, "\n")
## Dispersion parameter: 2.340484

Karena nilai dari Dispersion > 1, dimana 2.3404484 > 1 maka overdispersion mungkin terjadi sehingga model alternatif seperti negatif binomial regresi dapat digunakan.

plot(fitted(model_poisson), residuals(model_poisson, type = "deviance"),
 main = "Residual vs Fitted Values",
 xlab = "Nilai Prediksi", ylab = "Residual Deviance",
 pch = 16, col = "red")
abline(h = 0, lty = 2)

# Visualisasi hasil
data_toko$Prediksi <- predict(model_poisson, type = "response")
par(mfrow = c(1, 2))
boxplot(Prediksi ~ Promosi, data = data_toko,
 main = "Prediksi Pengunjung vs Jumlah Promosi",
 xlab = "Jumlah Promosi", ylab = "Prediksi Jumlah Pengunjung",
 col = "lightyellow")
points(as.factor(data_toko$Promosi), data_toko$Pengunjung, col = "blue", pch = 16, cex = 0.5)
plot(data_toko$UkuranToko, data_toko$Pengunjung,
 main = "Pengunjung: Observasi vs Prediksi",
 xlab = "Ukuran Toko (m²)", ylab = "Jumlah Pengunjung",
 pch = 16, col = "purple")
points(data_toko$UkuranToko, data_toko$Prediksi, col = "red", pch = 16)
legend("topleft", legend = c("Observasi", "Prediksi"),
 col = c("purple", "red"), pch = 16)

Interpretasi :

  • Intersep (\(\hat{\beta}_0\)): Estimasi intersep (2.800) menunjukkan log dari rata-rata jumlah pengunjung ketika Promosi = 0 dan UkuranToko = 0. Dalam skala asli, \(exp(2.800) \approx 16.44\) pengunjung. Namun, karena UkuranToko = 0 tidak realistis, intersep lebih berfungsi sebagai baseline.

  • Promosi (\(\hat{\beta}_1\)): Koefisien (0.105) berarti setiap tambahan satu kampanye promosi meningkatkan log(\(\mu\)) sebesar 0.105. Dalam skala asli, ini meningkatkan jumlah pengunjung sebesar \(exp(0.105) \approx 1.111\) kali (peningkatan 11.1%). P-value rendah (misalnya, \(2.7 \times 10^{-5}\)) menunjukkan efek signifikan (\(p < 0.05\)).

  • UkuranToko (\(\hat{\beta}_2\)): Koefisien (0.002) berarti kenaikan 1 m² pada ukuran toko meningkatkan log(\(\mu\)) sebesar 0.002, atau jumlah pengunjung meningkat sebesar \(exp(0.002) \approx 1.002\) kali (peningkatan 0.2% per m²). P-value rendah (misalnya, \(3 \times 10^{-5}\)) menunjukkan efek signifikan.

Inferensi GLM

Inferensi dalam Generalized Linear Model (GLM) adalah proses untuk membuat kesimpulan statistik tentang parameter model berdasarkan data yang diamati. GLM memperluas model regresi linear untuk menangani distribusi respons yang tidak normal, seperti distribusi Poisson, Binomial, atau Bernoulli, dengan menggunakan fungsi link untuk menghubungkan prediktor linier dengan rata-rata respons. Dalam inferensi GLM, kita fokus pada ekspektasi dan varians estimator, pengujian hipotesis, serta alat seperti uji Wald dan uji Likelihood Ratio untuk mengevaluasi signifikansi parameter.

Ekspektasi dan Varians dalam GLM

Ekspektasi Estimator

Ekpektasi estimator parameter \(\hat{\beta}\) menunjukkan apakah estimator tersebut unbiased (tidak bias). Dalam GLM, estimator Maximum Likelihood Estimation (MLE) untuk \(\hat{\beta}\) bersifat asymptotically unbiased, artinya ekspektasinya mendekati nilai sebenarnya \(\beta\) ketika ukuran sampel \(n\) besar.

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

Untuk distribusi dalam keluarga eksponensial, ekspektasi respons \(Y\) dihitung berdasarkan fungsi link dan prediktor linier:

\[ E[Y_i] = \mu_i = g^{-1}(x_i^{\top} \beta) \]

di mana \(g(\cdot)\) adalah fungsi link, \(x_i\) adalah vektor prediktor untuk observasi ke-\(i\), dan \(\beta\) adalah vektor parameter.

Contoh dalam Regresi Poisson

Untuk regresi Poisson, fungsi link-nya adalah logaritma, sehingga:

\[ \log(\mu_i) = x_i^{\top} \beta \quad \text{dan} \quad \mu_i = \exp(x_i^{\top} \beta) \]

Ekspektasi respons adalah:

\[ E[Y_i] = \mu_i = \exp(x_i^{\top} \beta) \]

Varians Estimator

Varians estimator \(\beta\) mengukur presisi estimasi. Dalam GLM, varians estimator MLE diaproksimasi menggunakan matriks informasi Fisher:

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

di mana :

  • \(X\) adalah matriks desain (matriks prediktor)

  • \(W\) adalah matriks bobot diagonal dengan elemen \[ w_{ii} = \frac{1}{Var(Y_i)} \left(\frac{\partial \mu_i}{\partial \theta_i}\right)^2, \quad \eta_i = x_i^T \beta \] adalah prediktor linier.

Untuk distribusi tertentu, varians respons \(Y_i\) bergantung pada rata-rata \(\mu_i\) dan parameter dispersi \(\phi\): \[ Var(Y_i) = \phi V(\mu_i) \]

Contoh Varians dalam Distribusi :

  • Poisson: \(V(\mu_i) = \mu_i\), sehingga \(Var(Y_i) = \mu_i\) (karena \(\phi = 1\)).

  • Binomial: \(V(\mu_i) = \mu_i(1 - \mu_i)\), sehingga \(Var(Y_i) = \frac{\mu_i(1 - \mu_i)}{n_i}\) untuk jumlah percobaan \(n_i\).

Perhitungan Varians Estimator: Matriks \(W\) untuk regresi logistik (binomial dengan link logit) memiliki elemen: \[ w_{ii} = \pi_i (1 - \pi_i) \] di mana \(\pi_i = \frac{e^{x_i^T \beta}}{1 + e^{x_i^T \beta}}\) adalah probabilitas sukses.

Contoh regresi poisson dalam syntax R :

# Simulasi data
set.seed(36)
x <- rnorm(100)
mu <- exp(0.5 + 0.8 * x)
y <- rpois(100, mu)
model <- glm(y ~ x, family = poisson)
summary(model)
## 
## Call:
## glm(formula = y ~ x, family = poisson)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.49463    0.08591   5.757 8.54e-09 ***
## x            0.81293    0.08015  10.143  < 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: 209.52  on 99  degrees of freedom
## Residual deviance: 101.83  on 98  degrees of freedom
## AIC: 321.69
## 
## Number of Fisher Scoring iterations: 5

Distribusi Asimptotik Estimator

Untuk ukuran sampel besar, estimator MLE \(\hat{\beta}\) mengikuti distribusi normal secara asimptotik:

\[ \hat{\beta} \sim \mathcal{N} \left( \beta, \left( X^\top W X \right)^{-1} \right) \]

Distribusi ini menjadi dasar untuk: - Pengujian hipotesis (misalnya, uji Wald). - Perhitungan p-value.

Interval Kepercayaan

Interval kepercayaan 95% untuk parameter \(\beta_j\) dihitung sebagai:

\[ \hat{\beta_j} \pm z_{0.025} \cdot SE(\hat{\beta_j}) \]

di mana: - \(z_{0.025} \approx 1.96\) adalah nilai kritis distribusi normal standar, - \(SE(\hat{\beta_j}) = \sqrt{\left( X^\top W X \right)^{-1}_{jj}}\) adalah standar error.

Contoh perhitungan interval kepercayaan :

# Simulasi data
set.seed(123)
n <- 100
x <- rnorm(n)
X <- cbind(1, x) # Matriks desain dengan intercept
beta_true <- c(0.5, 0.8)
eta <- X %*% beta_true
lambda <- exp(eta)
y <- rpois(n, lambda)
# Fit model Poisson
model <- glm(y ~ x, family = poisson)
model
## 
## Call:  glm(formula = y ~ x, family = poisson)
## 
## Coefficients:
## (Intercept)            x  
##      0.4495       0.8600  
## 
## Degrees of Freedom: 99 Total (i.e. Null);  98 Residual
## Null Deviance:       245 
## Residual Deviance: 106.8     AIC: 325.8
summary_model <- summary(model)
summary_model
## 
## Call:
## glm(formula = y ~ x, family = poisson)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  0.44950    0.08872   5.066 4.05e-07 ***
## x            0.86000    0.07463  11.523  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 245.05  on 99  degrees of freedom
## Residual deviance: 106.78  on 98  degrees of freedom
## AIC: 325.76
## 
## Number of Fisher Scoring iterations: 5
coef_table <- summary_model$coefficients

# Ekstrak estimasi dan standar error
beta_hat <- coef_table[, 1]
se_beta <- coef_table[, 2]
# Interval kepercayaan 95%
z <- qnorm(0.975)
ci_lower <- beta_hat - z * se_beta
ci_upper <- beta_hat + z * se_beta
# Tampilkan hasil
results <- data.frame(
 Parameter = names(beta_hat),
 Estimate = beta_hat,
 SE = se_beta,
 CI_Lower = ci_lower,
 CI_Upper = ci_upper
)
results

Interpretasi :

  • Nilai beta dengan nilai (0.44950 dan 0.8600) mendekati nilai sebenarnya

  • Nilai interval kepercayaan yang diperoleh pada tingkat kepercayaan 95% berada di rentang 0.2756061 < x < 0.6233841

# Visualisasi interval kepercayaan
library(ggplot2)
ggplot(results, aes(x = Parameter, y = Estimate)) +
 geom_point(size = 3) +
 geom_errorbar(aes(ymin = CI_Lower, ymax = CI_Upper), width = 0.2) +
 theme_minimal() +
 labs(title = "Interval Kepercayaan 95% untuk Parameter",
 y = "Estimasi", x = "Parameter") +
 coord_flip()

Interpretasi :

  • Plot menunjukkan estimasi parameter (titik) dan interval kepercayaan 95% (garis error).

  • Interval kepercayaan untuk x dan (Intercept) menunjukkan rentang nilai yang mungkin untuk parameter dengan probabilitas 95%.

  • Lebar interval mencerminkan presisi estimasi: interval yang lebih sempit menunjukkan standar error yang lebih kecil

Mencari Ekspektasi dan Varians dalam GLM

Ekspektasi dalam GLM dimodelkan melalui fungsi tautan yang menghubungkan rata-rata respons dengan kombinasi linier prediktor, sedangkan varians ditentukan oleh distribusi respons dan dinyatakan sebagai fungsi dari rata-rata, dengan parameter skala untuk fleksibilitas. Pemahaman yang tepat tentang struktur rata-rata dan varians sangat penting untuk memastikan model yang valid, terutama dengan penyesuaian untuk overdispersi atau data spesifik. Informasi ini didukung oleh sumber-sumber kredibel seperti buku teks dan kursus statistik universitas, memberikan dasar yang kuat untuk analisis data menggunakan GLM.

Ekspektasi

Ekspektasi, atau rata-rata kondisional \(E[Y|X=x]\), adalah fokus utama dalam GLM. Hubungan antara rata-rata dan prediktor diberikan oleh: \[ g(E[Y|X=x]) = \beta_0 + \beta_1x_1 + \cdots + \beta_p x_p \]Di sini, \(g\) adalah fungsi tautan, dan untuk menemukan ekspektasi, kita menerapkan invers fungsi tautan: \[ E[Y|X=x] = g^{-1}(\beta_0 + \beta_1x_1 + \cdots + \beta_p x_p) \]

Contoh spesifik:

Dalam regresi linier (Gaussian GLM), fungsi tautan adalah identitas, sehingga: \[ E[Y|X=x] = \beta_0 + \beta_1x_1 + \cdots + \beta_p x_p \]dengan distribusi normal.

Dalam regresi logistik (binomial GLM), fungsi tautan adalah logit, sehingga: \[ E[Y|X=x] = \frac{1}{1 + e^{-(\beta_0 + \beta_1x_1 + \cdots + \beta_p x_p)}} \]dengan distribusi binomial.

Dalam regresi Poisson, fungsi tautan adalah log, sehingga: \[ E[Y|X=x] = e^{\beta_0 + \beta_1x_1 + \cdots + \beta_p x_p} \]dengan distribusi Poisson.

Varians

Varians dalam GLM tidak diasumsikan homogen seperti dalam regresi linier biasa. Sebaliknya, varians ditentukan oleh distribusi variabel respons dan dinyatakan sebagai: \[ \text{Var}[Y|X=x] = \phi \cdot V(E[Y|X=x]) \] di mana:

\(\phi\) adalah parameter skala, yang sering kali diasumsikan 1 untuk distribusi standar, tetapi dapat diestimasi untuk overdispersi,

\(V(\mu)\) adalah fungsi varians, yang merupakan fungsi dari rata-rata \(\mu = E[Y|X=x]\).

Metode Penaksiran Parameter

Generalized Linear Models (GLM) digunakan untuk memodelkan data dengan distribusi seperti Poisson atau binomial. Untuk menyesuaikan model ini dengan data, kita perlu menaksir parameter, seperti koefisien yang menghubungkan variabel prediktor dengan hasil. Metode utama adalah Maximum Likelihood Estimation (MLE), yang mencari parameter yang membuat data paling mungkin terjadi. Karena solusi langsung sering tidak ada, metode iteratif seperti Newton-Raphson, Fisher scoring, dan IRLS digunakan untuk menemukan nilai parameter terbaik.

Maximum Likelihood Estimator (MLE)

MLE adalah metode utama untuk menaksir parameter dalam GLM. Tujuannya adalah memaksimalkan fungsi likelihood, yang mengukur probabilitas data yang diamati diberikan parameter tertentu. Dalam GLM, fungsi likelihood ditentukan oleh distribusi variabel respons, seperti distribusi Poisson untuk data hitungan atau binomial untuk data biner. Secara matematis, untuk observasi \(y_i\) dengan rata-rata yang diharapkan \(\mu_i\), fungsi log-likelihood \(\ell(\beta)\) didefinisikan berdasarkan distribusi yang dipilih. Misalnya, untuk regresi Poisson, log-likelihood adalah:

\[ \ell(\beta) = \sum_i \left[ y_i \log(\mu_i) - \mu_i - \log(y_i!) \right] \]

di mana \(\mu_i = e^{\mathbf{x}_i^\top \beta}\). Karena tidak ada solusi analitik untuk memaksimalkan \(\ell(\beta)\), metode iteratif diperlukan untuk menemukan \(\beta\) yang memaksimalkan likelihood.

Newton-Raphson

Metode Newton-Raphson adalah teknik iteratif untuk menemukan akar persamaan atau memaksimalkan fungsi, seperti log-likelihood dalam MLE. Dalam GLM, metode ini digunakan untuk menyelesaikan persamaan skor \(\ell'(\beta) = 0\), di mana \(\ell'(\beta)\) adalah gradien (turunan pertama) dari log-likelihood. Rumus pembaruan parameter adalah:

\[ \beta_t = \beta_{t-1} - \left[ \ell''(\beta_{t-1}) \right]^{-1} \ell'(\beta_{t-1}) \]

di mana \(\ell''(\beta)\) adalah matriks Hessian (turunan kedua). Proses ini dimulai dengan tebakan awal \(\beta_0\), sering kali nol, dan diperbarui hingga konvergensi, yaitu ketika perubahan parameter atau gradien cukup kecil (misalnya, kurang dari toleransi \(10^{-8}\)). Namun, metode ini bisa tidak stabil jika matriks Hessian tidak dapat diinversi atau jika data memiliki pola yang kompleks, seperti dalam kasus non-kanonik atau data dengan varians tinggi.

Fisher Scoring

Fisher scoring adalah variasi dari Newton-Raphson yang menggunakan matriks informasi Fisher (harapan dari turunan kedua negatif log-likelihood) sebagai ganti matriks Hessian. Rumus pembaruan adalah:

\[ \beta_t = \beta_{t-1} - \left[ \mathbb{E}[\ell''(\beta_{t-1})] \right]^{-1} \ell'(\beta_{t-1}) \]

Matriks informasi Fisher, \(\mathcal{I} = -\mathbb{E}[\ell''(\beta)]\), dijamin positif definit untuk distribusi dalam keluarga eksponensial kanonik, seperti Poisson atau binomial, sehingga lebih stabil daripada Hessian yang diamati. Dalam distribusi kanonik, Fisher scoring identik dengan Newton-Raphson karena \(\ell''(\beta) = \mathbb{E}[\ell''(\beta)]\). Keunggulan Fisher scoring adalah kemudahan perhitungan dan stabilitas numerik, terutama untuk data dengan distribusi standar.

IRLS (Iteratively Reweighted Least Square)

IRLS adalah algoritma spesifik untuk GLM yang setara dengan Fisher scoring dalam banyak kasus. IRLS mengubah masalah MLE menjadi serangkaian masalah weighted least squares yang diselesaikan secara iteratif. Rumus pembaruan dalam IRLS adalah:

\[ \beta_{t+1} = (\mathbf{X}^\top \mathbf{W}_t \mathbf{X})^{-1} \mathbf{X}^\top \mathbf{W}_t (\mathbf{G}_t (\mathbf{Y} - \mu_t) + \mathbf{X}\beta_t) \]

di mana \(\mathbf{X}\) adalah matriks prediktor, \(\mathbf{W}_t\) adalah matriks bobot berdasarkan varians \(\mu_t\), \(\mathbf{G}_t\) adalah matriks turunan fungsi tautan, \(\mathbf{Y}\) adalah vektor respons, dan \(\mu_t\) adalah rata-rata yang diprediksi pada iterasi ke-\(t\). IRLS efektif karena memanfaatkan struktur linier dari GLM, membuatnya mudah diimplementasikan dalam perangkat lunak seperti R atau SAS. Contohnya, dalam regresi logistik, bobot dihitung sebagai \(\mu_i (1 - \mu_i)\), dan proses berulang hingga konvergensi.

Diagnostik Model GLM

Analisis Residual

Salah satu metode diagnostik paling dasar dalam GLM adalah analisis residual, yang bertujuan untuk mengidentifikasi outlier, memeriksa kecocokan model, dan mendeteksi pola yang tidak dijelaskan oleh model. Residual Pearson, yang didefinisikan sebagai

\[ r_i = \frac{y_i - \hat{\mu}_i}{\sqrt{\hat{V}(Y_i)}} \]

digunakan untuk mendeteksi outlier dan menguji kecocokan model, di mana \(y_i\) adalah nilai observasi, \(\hat{\mu}_i\) adalah nilai yang diprediksi, dan \(\hat{V}(Y_i)\) adalah varians yang diperkirakan berdasarkan distribusi model. Dalam GLM, varians ini bergantung pada distribusi, misalnya, untuk Poisson,

\[ \hat{V}(Y_i) = \hat{\mu}_i \]

dan untuk binomial,

\[ \hat{V}(Y_i) = n_i \hat{\pi}_i (1 - \hat{\pi}_i) \]

Plot residual vs. prediktor linier, di mana prediktor linier adalah

\[ \hat{\eta}_i = g(\hat{\mu}_i) \]

dengan \(g\) sebagai fungsi tautan, dapat menunjukkan apakah ada pola yang tidak dijelaskan oleh model, seperti curvature, yang mungkin menunjukkan bahwa model rata-rata tidak sesuai. Misalnya, dalam regresi logistik, plot ini bisa menunjukkan apakah asumsi linearitas pada skala logit terpenuhi. Selain itu, deviance residual, yang merupakan generalisasi dari residual dalam model linier, digunakan untuk mendeteksi observasi yang tidak dijelaskan dengan baik oleh model. Deviance residual dihitung berdasarkan perbedaan antara log-likelihood observasi dan log-likelihood model penuh, dan sering digunakan dalam plot diagnostik untuk mengidentifikasi outlier.

# Hitung residual Pearson
residuals_pearson <- residuals(model, type = "pearson")
# Plot residual
ggplot(data.frame(Index = 1:n, Residual = residuals_pearson),
 aes(x = Index, y = Residual)) +
 geom_point(pch = 19) +
 geom_hline(yintercept = 0, color = "red", linetype = "dashed") +
 theme_minimal() +
 labs(title = "Plot Residual Pearson",
 x = "Indeks Observasi", y = "Residual Pearson")

Interpretasi :

  • Residual tersebar acak di sekitar garis nol tanpa pola tertentu, menunjukkan bahwa model Poisson sesuai dengan data.

  • Tidak ada pola sistematis seperti bentuk parabola, sehingga tidak ada indikasi pelanggaran asumsi model.

Pemeriksaan Kecocokan Model

Pemeriksaan dilakukan dengan membandingkan Plot Residual dengan Fitted Values

# Prediksi nilai fitted
fitted_values <- fitted(model)
# Plot residual vs fitted
ggplot(data.frame(Fitted = fitted_values, Residual = residuals_pearson),
 aes(x = Fitted, y = Residual)) +
 geom_point(pch = 19) +
 geom_hline(yintercept = 0, color = "red", linetype = "dashed") +
 theme_minimal() +
 labs(title = "Residual Pearson vs Nilai Fitted",
 x = "Nilai Fitted", y = "Residual Pearson")

Interpretasi :

  • Plot ini menunjukkan apakah residual memiliki varians konstan terhadap nilai fitted.

  • Penyebaran acak tanpa pola menunjukkan bahwa asumsi model Poisson terpenuhi.

Goodness of Fit

Uji goodness-of-fit penting untuk mengevaluasi seberapa baik model GLM cocok dengan data. Dalam GLM, deviance adalah statistik yang mengukur perbedaan antara likelihood model yang disesuaikan dan likelihood model saturasi (model dengan parameter untuk setiap observasi). Rumus deviance adalah

\[ D = 2 \sum_i \left[ l(y_i; y_i) - l(y_i; \hat{\mu}_i) \right] \]

di mana \(l\) adalah log-likelihood. Deviance rendah menunjukkan bahwa model cocok dengan baik, dan uji deviance dapat dibandingkan dengan distribusi chi-square untuk menentukan signifikansi. Selain itu, Pearson chi-square, yang didefinisikan sebagai

\[ X^2 = \sum_i \frac{(y_i - \hat{\mu}_i)^2}{\hat{V}(Y_i)} \]

digunakan untuk menguji apakah model sesuai dengan data, terutama untuk memeriksa overdispersi atau underdispersi.

Untuk regresi logistik, uji Hosmer-Lemeshow sering digunakan, yang membagi data ke dalam kelompok berdasarkan probabilitas yang diprediksi dan membandingkan jumlah observasi yang diamati dengan yang diharapkan. Jika p-value dari uji ini lebih besar dari 0,05, model dianggap cocok. Dalam praktiknya, uji ini berguna untuk data biner atau proporsi, tetapi untuk data hitungan seperti Poisson, deviance dan Pearson chi-square lebih umum digunakan.

Uji Asumsi Model

Uji asumsi model penting untuk memastikan bahwa model GLM memenuhi asumsi-asumsi dasarnya, seperti pemilihan fungsi tautan yang tepat dan distribusi yang sesuai. Salah satu asumsi penting adalah linearitas pada skala prediktor linier, yang dapat diuji dengan plot residual vs. prediktor atau dengan menambahkan transformasi prediktor (misalnya, \(\sqrt{x}\), \(\log x\), \(x^2\)) dan melihat apakah meningkatkan kecocokan model. Uji fungsi tautan dapat dilakukan dengan menambahkan kuadrat prediktor linier sebagai prediktor baru; jika signifikan, mungkin menunjukkan bahwa fungsi tautan, seperti identitas, logit, atau log, tidak sesuai. Misalnya, dalam regresi logistik, jika plot empirical logits vs. covariate menunjukkan non-linearitas, mungkin perlu transformasi prediktor.

Selain itu, plot residual atau absolute residual vs. prediktor linier dapat digunakan untuk mendeteksi heteroskedastisitas, yaitu varians yang tidak konstan. Dalam GLM, varians diharapkan sesuai dengan fungsi varians distribusi, seperti \(\mu\) untuk Poisson atau \(\mu(1-\mu)\) untuk binomial. Jika ada pola dalam plot ini, mungkin menunjukkan overdispersi, yang dapat diatasi dengan model seperti quasi-Poisson atau binomial negatif.

Metode Estimasi dan Inferensi Logistik

Regresi logistik digunakan untuk memodelkan probabilitas dari variabel respons biner (0/1) berdasarkan satu atau lebih variabel prediktor. Estimasi parameter dilakukan menggunakan Maximum Likelihood Estimation (MLE) karena model tidak linear dalam parameternya.

Estimasi dengan Newton-Raphson

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

Newton Raphson Iterasi Manual

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

Inferensi Parameter

Uji Wald

Uji Wald menguji hipotesis nol bahwa parameter tertentu sama dengan nol:

\[ H_0 : \beta_j = 0 \quad \text{vs} \quad H_1 : \beta_j \neq 0 \]

Statistik uji Wald adalah:

\[ Z = \frac{\hat{\beta_j}}{SE(\hat{\beta_j})} \]

di mana \(Z\) mengikuti distribusi normal standar \(N(0, 1)\) secara asimptotik. Kuadrat statistik Wald mengikuti distribusi chi-square dengan 1 derajat kebebasan:

\[ W = Z^2 = \left( \frac{\hat{\beta_j}}{SE(\hat{\beta_j})} \right)^2 \sim \chi_1^2 \]

P-value dihitung sebagai:

\[ p = P(\chi_1^2 > W) \]

Jika p-value < 0.05, kita menolak \(H_0\), menyimpulkan bahwa \(\beta_j\) signifikan.

Contoh perhitungan Uji Wald :

coef_val <- coef(model)[2] # Estimasi beta_1
se_val <- summary(model)$coefficients[2, 2] # Standar error
wald_z <- coef_val / se_val
wald_chisq <- wald_z^2
p_value <- pchisq(wald_chisq, df = 1, lower.tail = FALSE)
# Tampilkan hasil
cat("Statistik Z:", round(wald_z, 4), "\n")
## Statistik Z: 11.5229
cat("Statistik Chi-Square:", round(wald_chisq, 4), "\n")
## Statistik Chi-Square: 132.777
cat("P-value:", format(p_value, digits = 4), "\n")
## P-value: 1.012e-30

Interpretasi :

  • Statistik Z : Diperoleh nilai statistik Z sebesar 11.5229 sehinga dapat diketahui bahwa nilai Z sangat tinggi dan parameter yang diuji sangat berbeda dari 0 dan berkemungkinan besar signifikan

  • Statistik Chi-Squae : Diperoleh nilai Chi-Square sebesar 132.77 , dapat diketahui bahwa nilai p-value sangat besar sehingga dapat disimpulkan bahwa model memiliki kecocokan yang sangat baik dengan data

Uji Likelihood Ratio

Uji Likelihood Ratio (LRT) membandingkan dua model: model penuh (dengan semua prediktor) dan model tereduksi (tanpa prediktor tertentu). Statistik devians dihitung sebagai:

\[ D = -2 \left[ \ell(\hat{\beta}_{\text{reduced}}) - \ell(\hat{\beta}_{\text{full}}) \right] \]

di mana \(\ell(\cdot)\) adalah log-likelihood. Statistik \(D\) mengikuti distribusi chi-square dengan derajat kebebasan sama dengan selisih jumlah parameter antara kedua model.

Hipotesis :

  • \(H_0\): Model tereduksi cukup (prediktor tambahan tidak signifikan).

  • \(H_1\): Model penuh lebih baik.

# Model tereduksi (hanya intercept)
model_null <- glm(y ~ 1, family = poisson)
# Uji Likelihood Ratio
lrt <- anova(model_null, model, test = "Chisq")

#Tampilkan hasil
model_null
## 
## Call:  glm(formula = y ~ 1, family = poisson)
## 
## Coefficients:
## (Intercept)  
##      -1.079  
## 
## Degrees of Freedom: 99 Total (i.e. Null);  99 Residual
## Null Deviance:       73.36 
## Residual Deviance: 73.36     AIC: 143.4
lrt

Penerapan Regresi Logistik

set.seed(123)
n <- 100
x <- rnorm(n)
log_odds <- -0.5 + 1.2 * x
p <- 1 / (1 + exp(-log_odds))
y <- rbinom(n, 1, p)
data <- data.frame(x = x, y = y)
# Fit model logistik
model_logit <- glm(y ~ x, family = binomial)
# Ringkasan model
summary(model_logit)
## 
## Call:
## glm(formula = y ~ x, family = binomial)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -0.3097     0.2296  -1.349    0.177    
## x             1.2663     0.3080   4.111 3.94e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 137.99  on 99  degrees of freedom
## Residual deviance: 114.76  on 98  degrees of freedom
## AIC: 118.76
## 
## Number of Fisher Scoring iterations: 4

Uji Likelihood Ratio

# Model null
model_null_logit <- glm(y ~ 1, family = binomial, data = data)
# Uji Likelihood Ratio
anova(model_null_logit, model_logit, test = "Chisq")

Visualisasi Efek Prediktor

x_seq <- seq(min(x), max(x), length.out = 100)
new_data <- data.frame(x = x_seq)
pred_prob <- predict(model_logit, newdata = new_data, type = "response")
ggplot() +
 geom_point(aes(x = x, y = y), data = data, alpha = 0.5) +
 geom_line(aes(x = x_seq, y = pred_prob), color = "blue", size = 1) +
 theme_minimal() +
 labs(title = "Probabilitas Prediksi vs Prediktor",
 x = "x", y = "Probabilitas (y = 1)")
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

Keterangan Plot :

  • Titik-titik menunjukkan data observasi (0 atau 1)

  • Garis biru menunjukkan probabilitas prediksi dari model logistik.

  • Kurva sigmoid mencerminkan hubungan non-linear antara prediktor dan probabilitas sukses.

Keterangan Hasil :

  • Uji Wald menunjukkan bahwa prediktor x signifikan (p-value = 3.94e-05).

  • Uji Likelihood Ratio mengkonfirmasi bahwa model dengan prediktor x lebih baik daripada model null (pvalue = 1.438e-06).

  • Visualisasi probabilitas prediksi menunjukkan bagaimana memengaruhi peluang kejadian y = 1.

Evaluasi Kebaikan Model

Akaike Information Criterion (AIC)

AIC diperkenalkan oleh Hirotugu Akaike pada tahun 1974 sebagai alat untuk membandingkan model statistik berdasarkan teori informasi. AIC mengukur jumlah informasi yang hilang ketika model digunakan untuk merepresentasikan proses pembangkit data yang sebenarnya, dengan fokus pada meminimalkan divergensi Kullback-Leibler antara model dan kebenaran dasar. Rumus AIC didefinisikan sebagai:

\[ \text{AIC} = 2k - 2\ln(L) \]

di mana:

  • \(k\) adalah jumlah parameter dalam model (termasuk intersep dan parameter varians, jika ada)

  • \(L\) adalah nilai likelihood maksimum, yang mengukur seberapa baik model menjelaskan data yang diamati.

Komponen \(2k\) dalam rumus AIC adalah penalti untuk kompleksitas model, yang meningkat seiring bertambahnya jumlah parameter. Komponen \(-2\ln(L)\) mencerminkan kecocokan model, di mana nilai likelihood yang lebih tinggi (atau \(-2\ln(L)\) yang lebih rendah) menunjukkan kecocokan yang lebih baik. Model dengan nilai AIC terendah dianggap terbaik karena menyeimbangkan kecocokan dan kompleksitas.

AIC sangat berguna dalam konteks prediksi, seperti dalam peramalan deret waktu atau analisis regresi, karena mengestimasi kesalahan prediksi di luar sampel (out-of-sample prediction error). Namun, AIC cenderung memilih model yang lebih kompleks ketika jumlah parameter mendekati jumlah data, karena penalti \(2k\) relatif ringan.

aic <- AIC(model)
cat("AIC:", round(aic, 4), "\n")
## AIC: 325.7582

Bayesian Information Criterion (BIC)

Bayesian Information Criterion (BIC), yang diusulkan oleh Gideon Schwarz pada tahun 1978, bertujuan untuk memilih model terbaik dengan pendekatan yang berbeda, yaitu dengan fokus pada probabilitas posterior model dalam kerangka Bayesian. Rumus BIC adalah:

\[ \text{BIC} = k \ln(n) - 2\ln(L) \]

di mana:

  • \(k\) adalah jumlah parameter,

  • \(n\) adalah jumlah observasi (ukuran sampel),

  • \(L\) adalah likelihood maksimum.

Berbeda dengan AIC, BIC menggunakan penalti \(k \ln(n)\), yang bergantung pada ukuran sampel \(n\). Ketika \(n\) besar, penalti ini menjadi lebih berat dibandingkan penalti AIC (\(2k\)), sehingga BIC cenderung memilih model yang lebih sederhana.

BIC dirancang untuk mengaproksimasi probabilitas marginal data di bawah model, yang berguna untuk mengidentifikasi model yang benar jika model tersebut ada dalam himpunan kandidat. Karena sifatnya yang konsisten, BIC lebih disukai dalam situasi di mana tujuan utama adalah menemukan model yang paling mendekati kebenaran, terutama pada dataset besar.

bic <- BIC(model)
cat("BIC:", round(bic, 4), "\n")
## BIC: 330.9685

Visualisasi AIC dan BIC

# Data frame untuk visualisasi
criteria <- data.frame(
 Criterion = c("AIC", "BIC"),
 Value = c(aic, bic)
)
# Plot batang
ggplot(criteria, aes(x = Criterion, y = Value, fill = Criterion)) +
 geom_bar(stat = "identity") +
 theme_minimal() +
 labs(title = "Perbandingan AIC dan BIC",
 y = "Nilai", x = "Kriteria") +
 scale_fill_manual(values = c("AIC" = "blue", "BIC" = "red"))

Contoh Penerapan Regresi Logistik Multiple

library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
data(mtcars)
mtcars <- mtcars %>%
mutate(mpg_high = ifelse(mpg > median(mpg), 1, 0)) %>%
select(mpg_high, wt, hp, am) # pilih beberapa prediktor
head(mtcars)

Model Regresi Logistik Fit Model

model_logit <- glm(mpg_high ~ wt + hp + am, data = mtcars, family = binomial)
## Warning: glm.fit: algorithm did not converge
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: algorithm did not converge
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
summary(model_logit)
## 
## Call:
## glm(formula = mpg_high ~ wt + hp + am, family = binomial, data = mtcars)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)
## (Intercept)  6.244e+02  4.856e+05   0.001    0.999
## wt          -1.847e+02  1.875e+05  -0.001    0.999
## hp          -8.852e-02  2.045e+03   0.000    1.000
## am          -6.038e+01  1.044e+05  -0.001    1.000
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 4.4236e+01  on 31  degrees of freedom
## Residual deviance: 3.4224e-09  on 28  degrees of freedom
## AIC: 8
## 
## Number of Fisher Scoring iterations: 25

Uji Signifkansi Model (G² atau Likelihood Ratio Test)

model_null <- glm(mpg_high ~ 1, data = mtcars, family = binomial)
anova(model_null, model_logit, test = "Chisq")

Statistik Wald untuk Koefisien Individu

library(broom)
tidy(model_logit)

Interpretasi Odds Ratio Odds Ratio dan Confdence Interval

coefs <- summary(model_logit)$coefficients
OR <- exp(coefs[, "Estimate"])
CI_lower <- exp(coefs[, "Estimate"] - 1.96 * coefs[, "Std. Error"])
CI_upper <- exp(coefs[, "Estimate"] + 1.96 * coefs[, "Std. Error"])
odds_ratio_table <- cbind(OR, CI_lower, CI_upper)
round(odds_ratio_table, 3)
##                        OR CI_lower CI_upper
## (Intercept) 1.552814e+271        0      Inf
## wt           0.000000e+00        0      Inf
## hp           9.150000e-01        0      Inf
## am           0.000000e+00        0      Inf

Goodness of Fit Hosmer-Lemeshow Test

library(ResourceSelection)
## ResourceSelection 0.3-6   2023-06-27
hoslem.test(mtcars$mpg_high, fitted(model_logit), g=10)
## Warning in hoslem.test(mtcars$mpg_high, fitted(model_logit), g = 10): The data
## did not allow for the requested number of bins.
## 
##  Hosmer and Lemeshow goodness of fit (GOF) test
## 
## data:  mtcars$mpg_high, fitted(model_logit)
## X-squared = 1.7112e-09, df = 2, p-value = 1

Kesimpulan

Model logistik menunjukkan bahwa variabel wt dan am berpengaruh signifkan terhadap efsiensi bahan bakar (mpg_high). Odds ratio memberikan pemahaman intuitif terhadap kekuatan dan arah pengaruh prediktor.

Regresi Logistik dengan Prediktor Nominal, Ordinal, dan Rasio

Regresi logistik adalah metode statistik yang digunakan untuk memodelkan hubungan antara variabel dependen kategorikal, biasanya biner (misalnya, ya/tidak, lulus/gagal), dengan satu atau lebih variabel independen (prediktor). Prediktor dalam regresi logistik dapat berupa variabel nominal, ordinal, atau rasio, dan cara penanganan serta interpretasinya berbeda tergantung pada jenis prediktor tersebut.

Nominal

Prediktor nominal adalah variabel kategorikal yang tidak memiliki urutan alami, seperti jenis kelamin (laki-laki/perempuan), warna (merah/biru/hijau), atau golongan darah (A/B/AB/O). Dalam regresi logistik, prediktor nominal harus dikodekan menggunakan dummy coding (atau one-hot encoding) untuk dimasukkan ke dalam model. Proses ini melibatkan pembuatan variabel biner (0/1) untuk setiap kategori, kecuali satu kategori yang dijadikan referensi.

Contoh praktis dapat dilihat dalam analisis data, seperti memodelkan probabilitas seseorang untuk membeli produk berdasarkan “status pekerjaan” (karyawan/wiraswasta/pengangguran). Variabel ini dikodekan sebagai dummy, dan koefisien menunjukkan bagaimana status pekerjaan tertentu memengaruhi log-odds dibandingkan kategori referensi, seperti “pengangguran”.

Ordinal

Prediktor ordinal adalah variabel kategorikal yang memiliki urutan alami, tetapi jarak antar kategori tidak tentu sama, seperti tingkat pendidikan (SD/SMP/SMA/S1/S2) atau tingkat kepuasan (sangat tidak puas/tidak puas/netral/puas/sangat puas). Dalam regresi logistik, penanganan prediktor ordinal tergantung pada jenis model yang digunakan. Jika variabel dependen biner, prediktor ordinal dapat dimasukkan langsung dengan mengkodekan kategori sebagai angka berurutan (misalnya, 1, 2, 3, 4, 5), memperlakukannya seperti variabel kontinu, dengan asumsi bahwa efeknya linear pada skala log-odds.

Contoh praktis dapat dilihat dalam memodelkan probabilitas seseorang untuk mengajukan beasiswa, dengan prediktor ordinal seperti “tingkat pendidikan orang tua” (tidak tamat SD/tamat SD/tamat SMA/S1/S2). Dalam regresi logistik ordinal, prediktor ini dimasukkan dengan kode berurutan, dan koefisien menunjukkan bagaimana peningkatan tingkat pendidikan orang tua memengaruhi peluang.

Rasio

Prediktor rasio adalah variabel kontinu yang memiliki titik nol sejati dan skala yang konsisten, seperti usia, pendapatan, atau berat badan. Dalam regresi logistik, prediktor rasio dapat dimasukkan langsung ke dalam model tanpa perlu koding tambahan, karena sifatnya yang kontinu. Model regresi logistik mengasumsikan hubungan linear antara prediktor rasio dan log-odds dari variabel dependen. Namun, jika hubungan ini tidak linear, transformasi seperti logaritma, kuadrat, atau akar dapat digunakan untuk memodelkan hubungan yang lebih sesuai.

Contoh praktis dapat dilihat dalam memodelkan probabilitas seseorang untuk memiliki penyakit jantung, dengan prediktor rasio seperti “tekanan darah” atau “indeks massa tubuh (BMI)”. Prediktor ini dimasukkan langsung, dan koefisien menunjukkan bagaimana peningkatan tekanan darah atau BMI memengaruhi peluang penyakit.

Penentuan Pemilihan Model

Studi Kasus

Sebuah bank ingin memprediksi apakah seorang pelanggan akan berhasil mendapatkan persetujuan kredit berdasarkan jenis kelamin, tingkat pendidikan, dan pendapatan bulanan. Variabel dependen adalah keberhasilan pengajuan kredit (biner: 1 untuk disetujui, 0 untuk ditolak). Prediktor yang digunakan adalah :

  • Jenis Kelamin (nominal): Pria atau Wanita.

  • Tingkat Pendidikan (ordinal): SMP, SMA, Sarjana, atau Pascasarjana.

  • Pendapatan Bulanan (rasio): Dalam juta rupiah, sebagai variabel kontinu.

# Mengatur seed untuk reproduktibilitas
set.seed(123)

# Menentukan jumlah observasi
n <- 500

# Membuat variabel prediktor
jenis_kelamin <- sample(c("Pria", "Wanita"), n, replace = TRUE, prob = c(0.5, 0.5))
pendidikan <- sample(c("SMP", "SMA", "Sarjana", "Pascasarjana"), n, replace = TRUE, prob = c(0.2, 0.4, 0.3, 0.1))
pendapatan <- rnorm(n, mean = 75, sd = 20)

# Menghitung log-odds untuk probabilitas keberhasilan kredit
logit_p <- -1.5 + 0.4 * (jenis_kelamin == "Wanita") + 0.6 * as.numeric(factor(pendidikan, ordered = TRUE, levels = c("SMP", "SMA", "Sarjana", "Pascasarjana"))) + 0.02 * pendapatan

# Menghitung probabilitas
p <- 1 / (1 + exp(-logit_p))

# Menghasilkan variabel dependen (keberhasilan kredit)
persetujuan_kredit <- rbinom(n, 1, p)

# Membuat data frame menggunakan tibble
library(tibble)
data_kredit <- tibble(persetujuan_kredit, jenis_kelamin, pendidikan, pendapatan)

# Menampilkan 6 baris pertama dari data
head(data_kredit)

Dataset berisi 500 observasi dengan variabel persetujuan kredit, jenis kelamin, pendidikan dan pendapatan

Eksplorasi Data

data_kredit %>%
  dplyr::group_by(persetujuan_kredit) %>%
  dplyr::summarise(
    Jumlah = dplyr::n(),
    Rata2_pendapatan = mean(pendapatan)
)

Interpretasi :

  • Dari total 500 pelanggan, 103 pelanggan (20.6%) memiliki pengajuan kredit yang ditolak (persetujuan_kredit = 0), sedangkan 397 pelanggan (79.4%) memiliki pengajuan kredit yang disetujui (persetujuan_kredit = 1).

  • Pelanggan dengan pengajuan kredit yang ditolak memiliki rata-rata pendapatan bulanan sebesar 69.32 juta rupiah.

  • Pelanggan dengan pengajuan kredit yang disetujui memiliki rata-rata pendapatan bulanan sebesar 76.42 juta rupiah.

  • Perbedaan rata-rata pendapatan antara kelompok yang disetujui dan ditolak adalah sekitar 7.1 juta rupiah (76.42 - 69.32). Ini menunjukkan bahwa pelanggan dengan pendapatan lebih tinggi cenderung memiliki peluang lebih besar untuk mendapatkan persetujuan kredit.

Perlakuan Variabel Ordinal

Treat Sebagai Nominal

library(dplyr)
data_kredit_nominal <- data_kredit %>%
  mutate(
    pendidikan = factor(pendidikan, levels = c("SMP", "SMA", "Sarjana", "Pascasarjana"))
  )

model_nominal <- glm(persetujuan_kredit ~ jenis_kelamin + pendidikan + pendapatan, 
                     data = data_kredit_nominal, 
                     family = binomial(link = "logit"))

summary(model_nominal)
## 
## Call:
## glm(formula = persetujuan_kredit ~ jenis_kelamin + pendidikan + 
##     pendapatan, family = binomial(link = "logit"), data = data_kredit_nominal)
## 
## Coefficients:
##                         Estimate Std. Error z value Pr(>|z|)    
## (Intercept)            -1.163763   0.495043  -2.351 0.018732 *  
## jenis_kelaminWanita     0.769178   0.236836   3.248 0.001163 ** 
## pendidikanSMA           0.710698   0.274155   2.592 0.009533 ** 
## pendidikanSarjana       1.504334   0.331218   4.542 5.58e-06 ***
## pendidikanPascasarjana  1.922514   0.575824   3.339 0.000842 ***
## pendapatan              0.018348   0.005934   3.092 0.001988 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 508.61  on 499  degrees of freedom
## Residual deviance: 462.46  on 494  degrees of freedom
## AIC: 474.46
## 
## Number of Fisher Scoring iterations: 5

Interpretasi Koefisien :

  • Jenis Kelamin Wanita : Diperoleh nilai estimasi sebesar 0.769178 dimana variabel laki-laki adalah reference dapat diartikan bahwa wanita memiliki log odds yang lebih tinggi diterima dibandingkan laki laki, dimana peluang persetujuan kredit untuk wanita adalah 2.16 kali lebih tinggi dibandingkan laki-laki, dengan variabel lain konstan

  • Pendidikan SMA : Diperoleh nilai estimasi sebesar 0.710698 dimana variabel pendidikan SMP adalah reference dapat diartikan bahwa pelanggan dengan pendidikan SMA memiliki log odds persetujuan kredit lebih tinggi dibandingkan pendidikan SMP, dimana peluang persetujuan kredit untuk pelanggan dengan pendidikan SMA adalah 2.04 kali lebih tinggi dibandingkan SMP

  • Pendidikan Sarjana : Diperoleh nilai estimasi sebesar 1.504334 dimana variabel pendidikan SMP adalah reference dapat diartikan bahwa pelanggan dengan pendidikan sarjana memiliki log odds persetujuan kredit lebih tinggi dibandingkan pendidikan SMP, dimana peluang persetujuan kredit untuk pelanggan dengan pendidikan sarjana adalah 4,50 kali lebih tinggi dibandingkan SMP

  • Pendidikan Pascasarjana : Diperoleh nilai estimasi sebesar 1.922514 dimana variabel pendidikan SMP adalah reference dapat diartikan bahwa pelanggan dengan pendidikan pascasarjana memiliki log odds persetujuan kredit lebih tinggi dibandingkan pendidikan SMP, dimana peluang persetujuan kredit untuk pelanggan dengan pendidikan pascasarjana adalah 6,84 kali lebih tinggi dibandingkan SMP

  • Pendapatan : Diperoleh nilai estimasi sebesar 0.018348 sehingga hal ini menunjukkan perubahan log-odds untuk setiap kenaikan 1 juta rupiah dalam pendapatan bulanan. Artinya, setiap kenaikan 1 juta rupiah dalam pendapatan meningkatkan peluang persetujuan kredit sebesar 1.85%

Interpretasi Kecocokan Model :

  • Residual Deviance : Diperoleh nilai sebesar 462.46 pada derajat bebas 494, hal ini menunjukkan deviasi model dengan prediktor, yang lebih rendah dari null deviance (penurunan 508.61 - 462.46 = 46.15), menunjukkan bahwa prediktor (jenis_kelamin, pendidikan, pendapatan) menjelaskan sebagian variabilitas dalam persetujuan kredit.

Kesimpulan Gabungan :

  • Wanita memiliki peluang persetujuan kredit yang jauh lebih tinggi (2.16 kali) dibandingkan pria, yang mungkin mencerminkan pola dalam data simulasi atau perbedaan dalam profil kredit

  • Tingkat pendidikan yang lebih tinggi (SMA, Sarjana, Pascasarjana) secara signifikan meningkatkan peluang persetujuan kredit dibandingkan SMP. Efeknya meningkat secara monoton (SMA: 2.04x, Sarjana: 4.50x, Pascasarjana: 6.84x), meskipun pendidikan diperlakukan sebagai nominal.

  • Pendapatan bulanan memiliki efek positif kecil tetapi signifikan

Treat Sebagai Rasio

data_kredit_numeric <- data_kredit %>%
  mutate(
    pendidikan_numeric = case_when(
      pendidikan == "SMP" ~ 1,
      pendidikan == "SMA" ~ 2,
      pendidikan == "Sarjana" ~ 3,
      pendidikan == "Pascasarjana" ~ 4
    )
  )

# Menjalankan regresi logistik dengan pendidikan sebagai numerik
model_numeric <- glm(persetujuan_kredit ~ jenis_kelamin + pendidikan_numeric + pendapatan, 
                     data = data_kredit_numeric, 
                     family = binomial(link = "logit"))

# Menampilkan ringkasan model
summary(model_numeric)
## 
## Call:
## glm(formula = persetujuan_kredit ~ jenis_kelamin + pendidikan_numeric + 
##     pendapatan, family = binomial(link = "logit"), data = data_kredit_numeric)
## 
## Coefficients:
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)         -1.865313   0.551605  -3.382 0.000721 ***
## jenis_kelaminWanita  0.774024   0.236259   3.276 0.001052 ** 
## pendidikan_numeric   0.710172   0.141419   5.022 5.12e-07 ***
## pendapatan           0.018323   0.005924   3.093 0.001980 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 508.61  on 499  degrees of freedom
## Residual deviance: 462.71  on 496  degrees of freedom
## AIC: 470.71
## 
## Number of Fisher Scoring iterations: 4

Interpretasi :

  • Jenis Kelamin Wanita : Diperoleh nilai estimasi sebesar 0.774042 dimana variabel laki-laki adalah reference dapat diartikan bahwa wanita memiliki log odds yang lebih tinggi diterima dibandingkan laki laki, dimana peluang persetujuan kredit untuk wanita adalah 2.17 kali lebih tinggi dibandingkan laki-laki, dengan variabel lain konstan

  • Pendidikan Numerik : Karena pendidikan_numeric adalah variabel numerik (SMP = 1, SMA = 2, Sarjana = 3, Pascasarjana = 4), koefisien ini menunjukkan perubahan log-odds untuk setiap kenaikan satu tingkat pendidikan. Nilai 0.710172 berarti setiap peningkatan satu tingkat (misalnya, dari SMP ke SMA atau SMA ke Sarjana) meningkatkan log-odds sebesar 0.710172. Artinya, setiap kenaikan satu tingkat pendidikan meningkatkan peluang persetujuan kredit sebesar 2.03 kali, dengan variabel lain konstan.

  • Pendapatan : Koefisien ini menunjukkan perubahan log-odds untuk setiap kenaikan 1 juta rupiah dalam pendapatan bulanan. Nilai 0.018323 berarti setiap peningkatan 1 juta rupiah meningkatkan log-odds sebesar 0.018323. Artinya, setiap kenaikan 1 juta rupiah dalam pendapatan meningkatkan peluang persetujuan kredit sebesar 1.85%, dengan variabel lain konstan

Perbandingan Model

list(
AIC_Nominal = AIC(model_nominal),
AIC_Numeric = AIC(model_numeric)
)
## $AIC_Nominal
## [1] 474.462
## 
## $AIC_Numeric
## [1] 470.7121

Interpretasi :

Pendekatan numerik untuk pendidikan lebih sederhana dibandingkan model nominal sebelumnya (yang menggunakan tiga koefisien untuk SMA, Sarjana, dan Pascasarjana), dan AIC yang lebih rendah (470.71 vs. 474.46) menunjukkan bahwa model ini lebih efisien

Kesimpulan Keseluruhan

  • Wanita memiliki peluang persetujuan kredit yang jauh lebih tinggi (2.17 kali) dibandingkan pria, yang konsisten dengan model sebelumnya (2.16 kali). Ini mungkin mencerminkan pola dalam data simulasi atau perbedaan dalam profil kredit

  • Setiap kenaikan tingkat pendidikan (misalnya, dari SMP ke SMA, atau Sarjana ke Pascasarjana) meningkatkan peluang persetujuan kredit sebesar 2.03 kali. Pendekatan numerik ini mengasumsikan efek linier pada skala log-odds, yang masuk akal karena pendidikan dalam simulasi dihasilkan dengan koefisien 0.6 per tingkat.

  • Pendapatan memiliki efek positif kecil tetapi signifikan, dengan setiap kenaikan 1 juta rupiah meningkatkan peluang sebesar 1.85%. Ini konsisten dengan koefisien simulasi (0.02) dan model sebelumnya (0.018348), serta tabel ringkasan sebelumnya yang menunjukkan rata-rata pendapatan lebih tinggi untuk kelompok disetujui (76.42 juta vs. 69.32 juta rupiah).

Goodness of Fit

nullmod <- glm(persetujuan_kredit ~ 1, data = data_kredit, family = binomial)
r2_nominal <- 1 - (logLik(model_nominal)/logLik(nullmod))
r2_numeric <- 1 - (logLik(model_numeric)/logLik(nullmod))
list(
McFadden_R2_Nominal = r2_nominal,
McFadden_R2_Numeric = r2_numeric
)
## $McFadden_R2_Nominal
## 'log Lik.' 0.09073084 (df=6)
## 
## $McFadden_R2_Numeric
## 'log Lik.' 0.09023921 (df=4)

Interpretasi :

  • Nilai McFadden R² sebesar 0.0907 menunjukkan bahwa model nominal menjelaskan sekitar 9.07% dari variabilitas dalam data persetujuan_kredit dibandingkan model nol.

  • Nilai McFadden R² sebesar 0.0902 menunjukkan bahwa model numerik menjelaskan sekitar 9.02% dari variabilitas dalam data persetujuan_kredit dibandingkan model nol.

  • Perbedaan McFadden R² antara kedua model sangat kecil (0.0907 - 0.0902 = 0.0005 atau 0.05%). Ini menunjukkan bahwa model numerik, meskipun lebih sederhana, memiliki kemampuan penjelasan yang hampir setara dengan model nominal. Dengan kata lain, asumsi efek linier pendidikan (dalam model numerik) tidak mengurangi kecocokan model secara signifikan dibandingkan dengan pendekatan nominal yang lebih fleksibel.

Visualisasi Prediksi

library(ggplot2)
# Menambahkan probabilitas prediksi ke data frame
data_kredit_nominal <- data_kredit_nominal %>% 
  mutate(predicted = predict(model_nominal, type = "response"))

data_kredit_numeric <- data_kredit_numeric %>% 
  mutate(predicted = predict(model_numeric, type = "response"))

# Memeriksa apakah kolom 'predicted' ada
str(data_kredit_nominal)
## tibble [500 × 5] (S3: tbl_df/tbl/data.frame)
##  $ persetujuan_kredit: int [1:500] 1 1 1 1 1 1 1 1 1 1 ...
##  $ jenis_kelamin     : chr [1:500] "Wanita" "Pria" "Wanita" "Pria" ...
##  $ pendidikan        : Factor w/ 4 levels "SMP","SMA","Sarjana",..: 2 2 2 2 2 2 3 3 4 2 ...
##  $ pendapatan        : num [1:500] 63 55.1 95.5 90 44.8 ...
##  $ predicted         : Named num [1:500] 0.813 0.636 0.888 0.768 0.591 ...
##   ..- attr(*, "names")= chr [1:500] "1" "2" "3" "4" ...
str(data_kredit_numeric)
## tibble [500 × 6] (S3: tbl_df/tbl/data.frame)
##  $ persetujuan_kredit: int [1:500] 1 1 1 1 1 1 1 1 1 1 ...
##  $ jenis_kelamin     : chr [1:500] "Wanita" "Pria" "Wanita" "Pria" ...
##  $ pendidikan        : chr [1:500] "SMA" "SMA" "SMA" "SMA" ...
##  $ pendapatan        : num [1:500] 63 55.1 95.5 90 44.8 ...
##  $ pendidikan_numeric: num [1:500] 2 2 2 2 2 2 3 3 4 2 ...
##  $ predicted         : Named num [1:500] 0.815 0.638 0.889 0.769 0.593 ...
##   ..- attr(*, "names")= chr [1:500] "1" "2" "3" "4" ...
# Plot untuk model nominal
plot_nominal <- data_kredit_nominal %>%
  ggplot(aes(x = pendapatan, y = predicted, color = pendidikan)) +
  geom_point(alpha = 0.6) +
  labs(title = "Prediksi Probabilitas (Pendidikan sebagai Nominal)", 
       x = "Pendapatan (Juta Rupiah)", 
       y = "Probabilitas Prediksi Persetujuan Kredit") +
  theme_minimal()

# Plot untuk model numerik
plot_numeric <- data_kredit_numeric %>%
  ggplot(aes(x = pendapatan, y = predicted, color = pendidikan)) +
  geom_point(alpha = 0.6) +
  labs(title = "Prediksi Probabilitas (Pendidikan sebagai Numerik)", 
       x = "Pendapatan (Juta Rupiah)", 
       y = "Probabilitas Prediksi Persetujuan Kredit") +
  theme_minimal()

# Menampilkan plot
print(plot_nominal)

print(plot_numeric)

Pemilihan Model Regresi Logistik dan Evaluasi

Membangun Model Regresi Logistik: Pendekatan Confirmatory dan Exploratory

Pendekatan Confirmatory

Pendekatan confirmatory (konfirmatori) adalah pendekatan berbasis hipotesis yang bertujuan untuk menguji teori atau hubungan spesifik yang telah ditentukan sebelumnya berdasarkan literatur atau penelitian sebelumnya. Dalam konteks regresi logistik, pendekatan ini digunakan untuk mengkonfirmasi apakah variabel independen tertentu memiliki efek signifikan terhadap variabel dependen biner, sesuai dengan hipotesis yang telah dirumuskan.

Karakteristik

  • Berbasis teori atau penelitian sebelumnya.
  • Variabel independen dipilih berdasarkan hipotesis yang jelas.
  • Fokus pada pengujian statistik untuk mengkonfirmasi atau menolak hipotesis.
  • Model spesifikasi (misalnya, variabel yang dimasukkan) ditentukan sebelum analisis data dilakukan.

Pendekatan Exploratory

Pendekatan exploratory (eksploratori) adalah pendekatan berbasis data yang bertujuan untuk menemukan pola atau hubungan yang tidak diketahui sebelumnya antara variabel independen dan variabel dependen. Pendekatan ini sering digunakan ketika teori atau hipotesis belum jelas, atau saat peneliti ingin mengidentifikasi prediktor potensial.

Karakteristik

  • Tidak bergantung pada hipotesis spesifik, tetapi lebih pada eksplorasi data.
  • Variabel independen dipilih berdasarkan analisis data awal, seperti korelasi atau uji univariat.
  • Fokus pada penemuan pola baru atau hubungan tak terduga.
  • Model sering dimodifikasi selama analisis (misalnya, dengan menambah/menghapus variabel).

Contoh Studi Kasus

Sebuah perusahaan ritel di Indonesia ingin meluncurkan produk baru (misalnya, minuman kesehatan berbasis herbal). Untuk mengoptimalkan anggaran pemasaran, perusahaan perlu memprediksi apakah seorang pelanggan potensial akan membeli produk tersebut berdasarkan data demografi dan perilaku mereka. Perusahaan telah mengumpulkan data dari 200 pelanggan potensial selama uji coba kampanye pemasaran awal.

Variabel - Variabel :

  • Pembelian (Y) : 1 = Membeli, 0 = Tidak

  • Usia (X1) : Usia pelanggan

  • Status Langganan (X2) : 1 = Pelanggan setia, 0 = Pelanggan baru

  • Frekuensi Kunjungan (X3) : Frekuensi kunjungan pelanggan ke toko atau situs web perusahaan dalam sebulan terakhir

# Mengatur seed untuk reproduktibilitas
set.seed(123)

# Menentukan jumlah observasi
n <- 200

# Membuat variabel independen
usia <- rnorm(n)  # Usia pelanggan (kontinu)
status_langganan <- rbinom(n, 1, 0.5)  # Status langganan (1 = Pelanggan Setia, 0 = Pelanggan Baru)
frekuensi_kunjungan <- rnorm(n)  # Frekuensi kunjungan per bulan (kontinu)

# Menghitung linier prediktor untuk probabilitas pembelian
lin_pred <- -0.5 + 1.2 * usia - 0.8 * status_langganan + 0.5 * frekuensi_kunjungan

# Menghitung probabilitas pembelian
p <- 1 / (1 + exp(-lin_pred))

# Menghasilkan variabel dependen (pembelian)
pembelian <- rbinom(n, 1, p)

# Membuat data frame
df_pemasaran <- data.frame(pembelian = as.factor(pembelian), usia, status_langganan, frekuensi_kunjungan)

# Menampilkan 6 baris pertama data frame
head(df_pemasaran)

Pemilihan Model

Model Full

model_full <- glm(pembelian ~ usia + status_langganan + frekuensi_kunjungan, data = df_pemasaran, family = binomial)
summary(model_full)
## 
## Call:
## glm(formula = pembelian ~ usia + status_langganan + frekuensi_kunjungan, 
##     family = binomial, data = df_pemasaran)
## 
## Coefficients:
##                     Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          -0.7148     0.2470  -2.894  0.00381 ** 
## usia                  1.4029     0.2315   6.061 1.35e-09 ***
## status_langganan     -0.2507     0.3463  -0.724  0.46903    
## frekuensi_kunjungan   0.3567     0.1704   2.094  0.03630 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 257.72  on 199  degrees of freedom
## Residual deviance: 202.67  on 196  degrees of freedom
## AIC: 210.67
## 
## Number of Fisher Scoring iterations: 4

Interpretasi :

  • Setiap kenaikan satu standar deviasi usia meningkatkan peluang pembelian sebesar 4.07 kali, dengan variabel lain (status_langganan dan frekuensi_kunjungan) konstan. P-value sangat kecil (< 0.001) menunjukkan efek yang sangat signifikan.

  • Peluang pembelian untuk pelanggan setia adalah 77.8% dari peluang pelanggan baru, dengan variabel lain konstan. Namun, p-value 0.469 (> 0.05) menunjukkan bahwa efek ini tidak signifikan secara statistik.

Metode Stepwise : Forward, Backward, Kedua Arah

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

AIC(model_full, step_forward, step_backward, step_both)

Interpretasi :

  • Ketiga model seleksi stepwise (step forward, step backward, step both) menghasilkan model yang sama, dengan 3 parameter dan AIC = 209.1998.

  • Jumlah parameter yang berbeda antara 3 dengan 4 menandakan bahwa satu variabel x telah dihapus dari model, dan variabel tersebut adalah status langganan dikarenakan berdasarkan informasi sebelumnya pada model full bahwa variabel tersebut tidak signifikan

  • Nilai AIC yang lebih kecil pada metode stepwise dibandingkan dengan nilai AIC pada model full menandakan bahwa metode stepwise lebih baik dibandingkan dengan model full, hal ini mungkin terjadi dikarenakan pada metode stepwise dilakukan penghapusan variabel yang tidak signifikan sehingga meningkatkan kekuatan model tersebut

Evaluasi Model : ROC dan AUC

library(pROC)
## Type 'citation("pROC")' for a citation.
## 
## Attaching package: 'pROC'
## The following objects are masked from 'package:stats':
## 
##     cov, smooth, var
pred_prob <- predict(step_both, type = "response")
roc_obj <- roc(df_pemasaran$pembelian, pred_prob)
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
## Setting levels: control = 0, case = 1
## Setting direction: controls < cases
plot(roc_obj, main = "Kurva ROC", col = "blue")

Interpretasi :

  • Kurva ROC yang berada di atas garis diagonal menunjukkan bahwa model regresi logistik (pembelian ~ usia + frekuensi_kunjungan) efektif dalam membedakan pelanggan yang membeli dari yang tidak membeli.
  • AUC yang diperkirakan 0.7–0.8 menunjukkan kemampuan diskriminasi yang baik. Dalam konteks pemasaran, ini berarti model dapat membantu mengidentifikasi pelanggan dengan probabilitas pembelian tinggi dengan akurasi yang memadai.

Pseudo R-Squared

McFadden R-Squared

McFadden R-squared, yang diperkenalkan oleh Daniel McFadden, mengukur perbaikan model logistik dibandingkan dengan model nol (hanya intersep) berdasarkan log-likelihood. Rumusnya adalah:

\[ R^2_{\text{McFadden}} = 1 - \frac{\ln(L_{\text{model}})}{\ln(L_{\text{null}})} \]

di mana:

  • \(\ln(L_{\text{model}})\) adalah log-likelihood dari model yang diusulkan (dengan prediktor),

  • \(\ln(L_{\text{null}})\) adalah log-likelihood dari model nol (tanpa prediktor, hanya intersep).

Interpretasi :

  • Nilai berkisar dari 0 hingga 1 (teoretis), tetapi biasanya lebih kecil daripada R-squared tradisional (0.2–0.4 sering dianggap baik untuk model logistik).
  • Nilai mendekati 0 menunjukkan bahwa model tidak lebih baik dari model nol, sedangkan nilai lebih tinggi menunjukkan model yang lebih baik menjelaskan data.
  • Contoh: Dalam studi kasus sebelumnya, McFadden R² untuk model nominal adalah 0.0907 dan untuk model numerik 0.0902, menunjukkan bahwa model menjelaskan sekitar 9% variabilitas dibandingkan model nol.

Cox & Snell R-Squared

Cox & Snell R-squared, yang dikembangkan oleh David Cox dan E. Jean Snell, adalah ukuran pseudo R-squared yang berdasarkan pada kemungkinan (likelihood) dan memperhitungkan distribusi binomial dari variabel dependen. Rumusnya adalah:

\[ R^2_{\text{Cox \& Snell}} = 1 - \left(\frac{L_{\text{null}}}{L_{\text{model}}}\right)^{2/n} \]

di mana:

  • \(L_{\text{null}}\) adalah likelihood dari model nol,

  • \(L_{\text{model}}\) adalah likelihood dari model yang diusulkan,

  • \(n\) adalah jumlah observasi.

Interpretasi :

  • Nilai berkisar dari 0 hingga nilai maksimum yang kurang dari 1 (tergantung pada sample size dan distribusi data), biasanya mendekati 0.7–0.8 sebagai batas atas.
  • Nilai lebih tinggi menunjukkan model yang lebih baik, tetapi nilai maksimum tidak selalu 1, yang membatasi interpretasinya.
  • Contoh: Dalam model regresi logistik, nilai ini sering digunakan sebagai indikator awal kecocokan, tetapi perlu dibandingkan dengan varian lain untuk konteks yang lebih jelas

Nagelkerke R-Squared

Nagelkerke R-squared adalah modifikasi dari Cox & Snell R-squared, yang dikembangkan oleh Jan Nagelkerke untuk menormalkan nilai sehingga batas maksimumnya adalah 1. Ini memungkinkan interpretasi yang lebih mirip dengan R-squared linier. Rumusnya adalah:

\[ R^2_{\text{Nagelkerke}} = \frac{1 - \left(\frac{L_{\text{null}}}{L_{\text{model}}}\right)^{2/n}}{1 - \left(L_{\text{null}}\right)^{2/n}} \]

Penyebut : \(1 - (L_{\text{null}})^{2/n}\) menormalkan nilai agar maksimumnya menjadi 1 ketika model sempurna.

Interpretasi :

  • Nilai berkisar dari 0 hingga 1, dengan 1 menunjukkan model yang sempurna dan 0 menunjukkan model nol.
  • Nilai 0.2–0.4 sering dianggap baik untuk model logistik, mirip dengan R-squared linier.
  • Contoh: Jika Nagelkerke R² = 0.15, model menjelaskan 15% variabilitas maksimum yang mungkin, yang lebih mudah dibandingkan dengan McFadden atau Cox & Snell.
library(DescTools)
PseudoR2(step_both, which = c("CoxSnell", "Nagelkerke", "McFadden"))
##   CoxSnell Nagelkerke   McFadden 
##  0.2385981  0.3294000  0.2115439

Interpretasi :

  • McFadden (0.2115) < Cox & Snell (0.2386) < Nagelkerke (0.3294), yang konsisten dengan sifat masing-masing varian. Nagelkerke selalu memberikan nilai tertinggi karena normalisasinya, sementara McFadden cenderung lebih konservatif.

  • Nilai pseudo R-squared yang berkisar antara 0.2115 (McFadden) dan 0.3294 (Nagelkerke) menunjukkan bahwa model menjelaskan 21.15% hingga 32.94% variabilitas dalam keputusan pembelian

Tabel Klasifikasi dan Evaluasi

library(caret)
## Loading required package: lattice
## 
## Attaching package: 'caret'
## The following objects are masked from 'package:DescTools':
## 
##     MAE, RMSE
pred_class <- ifelse(pred_prob >= 0.5, 1, 0)
conf_matrix <- confusionMatrix(factor(pred_class), df_pemasaran$pembelian, positive = "1")
conf_matrix
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction   0   1
##          0 116  32
##          1  15  37
##                                        
##                Accuracy : 0.765        
##                  95% CI : (0.7, 0.8219)
##     No Information Rate : 0.655        
##     P-Value [Acc > NIR] : 0.0005028    
##                                        
##                   Kappa : 0.4478       
##                                        
##  Mcnemar's Test P-Value : 0.0196041    
##                                        
##             Sensitivity : 0.5362       
##             Specificity : 0.8855       
##          Pos Pred Value : 0.7115       
##          Neg Pred Value : 0.7838       
##              Prevalence : 0.3450       
##          Detection Rate : 0.1850       
##    Detection Prevalence : 0.2600       
##       Balanced Accuracy : 0.7109       
##                                        
##        'Positive' Class : 1            
## 

Interpretasi :

  • Akurasi sebesar 0.765 (76.5%) menunjukkan bahwa 76.5% dari semua prediksi (baik 0 maupun 1) adalah benar. Interval kepercayaan 95% (0.7 hingga 0.8219) menunjukkan bahwa akurasi tinggi kemungkinan besar berada dalam rentang ini dengan tingkat kepercayaan 95%.

  • NIR adalah proporsi kelas mayoritas (dalam hal ini, 131 dari 200 pelanggan tidak membeli, atau 0.655). P-value < 0.05 (0.0005028) menunjukkan bahwa akurasi model (0.765) secara signifikan lebih baik daripada hanya menebak kelas mayoritas, membuktikan bahwa model memiliki nilai prediktif

  • Kappa mengukur kesepakatan di luar kesempatan acak. Nilai 0.4478 menunjukkan kesepakatan moderat (berdasar skala Landis & Koch: 0.41–0.60 = moderat). Ini mengindikasikan bahwa model lebih baik daripada tebakan acak, tetapi masih ada ruang untuk perbaikan.

  • Uji McNemar membandingkan kesalahan prediksi antara kelas. P-value < 0.05 (0.0196) menunjukkan adanya perbedaan signifikan dalam pola kesalahan

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

Interpretasi :

  • Sensitivitas 0.5362 (53.62%) adalah proporsi pelanggan yang benar-benar membeli (1) yang diklasifikasikan sebagai membeli oleh model.

  • Spesifisitas 0.8855 (88.55%) adalah proporsi pelanggan yang tidak membeli (0) yang diklasifikasikan sebagai tidak membeli.

Metode Perbandingan Model dalam Regresi Logistik

Pembuatan Model

model1 <- glm(pembelian ~ usia, data = df_pemasaran, family = binomial)
model2 <- glm(pembelian ~ usia + status_langganan, data = df_pemasaran, family = binomial)
model3 <- glm(pembelian ~ usia + status_langganan + frekuensi_kunjungan, data = df_pemasaran, family = binomial)

Perbandingan AIC dan Deviance

model_comp <- data.frame(
  Model = c("Model 1", "Model 2", "Model 3"),
  AIC = c(AIC(model1), AIC(model2), AIC(model3)),
  Deviance = c(deviance(model1), deviance(model2), deviance(model3))
)
model_comp

Interpretasi :

  • Model 3 memiliki AIC dan Deviance terendah di antara ketiga model, menunjukkan bahwa ini adalah model terbaik berdasarkan kedua metrik. Nilai Deviance 202.6730 lebih rendah dibandingkan null deviance (257.72 dari output sebelumnya), menunjukkan peningkatan signifikan dalam kecocokan. Karena AIC Model 3 (210.6739) mendekati AIC dari step_both sebelumnya (209.1998)

  • Model 1 memiliki AIC dan Deviance yang lebih tinggi daripada Model 3. Perbedaan AIC dengan Model 3 adalah 211.8410 - 210.6739 = 1.1671, yang < 2, menunjukkan bahwa Model 1 hampir setara dengan Model 3

  • Model 2 memiliki AIC dan Deviance tertinggi. Perbedaan AIC dengan Model 3 adalah 213.2377 - 210.6739 = 2.5638, yang berada di rentang 2–4, menunjukkan bahwa Model 3 sedikit lebih baik.

Likelihood Ratio-Test

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

Interpretasi :

  • Penurunan Deviance sebesar 0.6033 setelah menambahkan status langganan. Penurunan ini kecil, menunjukkan bahwa status langganan hanya memberikan kontribusi minimal terhadap penjelasan variabilitas dalam pembelian

  • P-value 0.4373 > 0.05 menunjukkan bahwa penambahan status langganan ke Model 1 untuk membentuk Model 2 tidak signifikan secara statistik. Ini berarti tidak ada bukti yang cukup untuk menyatakan bahwa status langganan meningkatkan kecocokan model secara bermakna dibandingkan hanya menggunakan usia

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

Intepretasi :

  • Penurunan Deviance sebesar 4.5638 setelah menambahkan frekuensi kunjungan ke Model. Penurunan ini menunjukkan bahwa frekuensi kunjungan membantu menjelaskan variabilitas tambahan dalam pembelian

  • P-value 0.03265 < 0.05 menunjukkan bahwa penambahan frekuensi kunjungan ke Model 1 untuk membentuk Model 2 signifikan secara statistik. Ini berarti frekuensi kunjungan memberikan kontribusi bermakna terhadap peningkatan kecocokan model.

Prinsip Parsimony

Prinsip parsimony menekankan bahwa model atau teori yang menggunakan jumlah asumsi, parameter, atau kompleksitas minimum, namun tetap mampu menjelaskan data dengan memadai. Ini didasarkan pada :

  • Model yang lebih sederhana cenderung lebih mudah diinterpretasikan dan diuji.
  • Kompleksitas tambahan harus dibenarkan oleh bukti yang jelas, bukan spekulasi.
  • Model yang lebih kompleks berisiko overfitting (menyesuaikan terlalu baik dengan data pelatihan tetapi gagal menggeneralisasi).

Penerapan Prinsip Parsimony

  • AIC (Akaike Information Criterion) : Karenamenentukan model berdasarkan jumlah parameter (kompleksitas) sambil mempertimbangkan kecocokan terhadap data.

  • Deviance : Karena mengukur seberapa jauh model saat ini dibandingkan dengan model sempurna

  • Likelihood Ratio-Test : Karena metode ini menguji apakah penambahan variabel dalam model secara signifikan meningkatkan kecocokan model

Pemilihan Treshold Optimal

# Mengatur seed untuk reproduktibilitas
set.seed(123)

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

# Menghitung sensitivitas untuk setiap threshold
results$Sensitivity <- sapply(thresholds, function(t) {
  pred_class <- ifelse(pred_prob >= t, 1, 0)  # Gunakan pred_prob
  cm <- table(Pred = pred_class, Obs = df_pemasaran$pembelian)  # Gunakan df_pemasaran$pembelian
  TP <- ifelse("1" %in% rownames(cm) && "1" %in% colnames(cm), cm["1", "1"], 0)  # Pengecekan TP
  FN <- ifelse("0" %in% rownames(cm) && "1" %in% colnames(cm), cm["0", "1"], 0)  # Pengecekan FN
  if (TP + FN > 0) TP / (TP + FN) else NA  # Hindari pembagian oleh nol
})

# Menampilkan hasil
print(results)
##    Threshold Sensitivity
## 1       0.10  0.98550725
## 2       0.15  0.92753623
## 3       0.20  0.84057971
## 4       0.25  0.81159420
## 5       0.30  0.75362319
## 6       0.35  0.69565217
## 7       0.40  0.68115942
## 8       0.45  0.62318841
## 9       0.50  0.53623188
## 10      0.55  0.47826087
## 11      0.60  0.42028986
## 12      0.65  0.33333333
## 13      0.70  0.30434783
## 14      0.75  0.18840580
## 15      0.80  0.13043478
## 16      0.85  0.05797101
## 17      0.90  0.02898551

Interpretasi :

  • Karena cut off optimal bisa dipilih berdasarkan Maksimum Senitivity + Specificity maka Threshold optimal diperkirakan di kisaran 0.35–0.40, di mana sensitivitas (0.6957–0.6812) dan spesifisitas (estimasi ~0.70–0.75) memberikan kombinasi terbaik.

  • Jika False Negatives Harus Dihindari maka kita akan memilih threshold 0.15 (sensitivitas 0.9275) untuk memaksimalkan deteksi pembeli, meskipun spesifisitas akan rendah.

  • Jika Efisiensi Diutamakan maka kita akan pilih threshold 0.50 (sensitivitas 0.5362), yang konsisten dengan matriks kebingungan sebelumnya dan memberikan spesifisitas tinggi (0.8855), mengoptimalkan targeting dengan mengurangi false positives.

Precision Recall Curve (PR Curve)

Kurva Precision-Recall adalah alat visual dan analitis yang digunakan dalam evaluasi model klasifikasi biner, seperti regresi logistik yang Anda gunakan untuk memprediksi keberhasilan pembelian produk kesehatan. Kurva ini menggambarkan hubungan antara precision (ketepatan) dan recall (panggilan kembali atau sensitivitas) pada berbagai nilai ambang batas (threshold) probabilitas prediksi.

  • Precision (Ketetapan)

Precision diukur sebagai rasio prediksi positif yang benar (True Positives, TP) terhadap total prediksi positif (TP + False Positives, FP). Rumusnya adalah:

\[ \text{Precision} = \frac{\text{TP}}{\text{TP} + \text{FP}} \]

Artinya: Dari semua instance yang diprediksi sebagai kelas positif (misalnya, “membeli”), berapa persen yang benar-benar positif. Precision mencerminkan kualitas prediksi positif.

  • Recall (Sensitivitas)

Recall diukur sebagai rasio prediksi positif yang benar (TP) terhadap total instance positif sebenarnya (TP + False Negatives, FN). Rumusnya adalah:

\[ \text{Recall} = \frac{\text{TP}}{\text{TP} + \text{FN}} \]

Artinya: Dari semua instance yang benar-benar positif, berapa persen yang berhasil dideteksi oleh model. Recall mencerminkan kemampuan model untuk menemukan semua kasus positif.

  • Trade off Precision dan Recall

Ketika threshold probabilitas digeser (misalnya, dari 0.9 ke 0.1), ada tradeoff: meningkatkan recall biasanya menurunkan precision, dan sebaliknya. Ini karena threshold yang lebih rendah akan mengklasifikasikan lebih banyak instance sebagai positif, meningkatkan TP dan FN, tetapi juga meningkatkan FP.

Contoh Kasus (Menggunakan Studi Kasus yang Sebelumnya) :

library(PRROC)
## Loading required package: rlang
# Mengatur seed untuk reproduktibilitas
set.seed(123)

# Membuat data simulasi sesuai kasus (mengganti x1, x2, x3 dengan variabel kasus)
usia <- rnorm(200)                    # Mengganti x1 dengan usia
status_langganan <- rbinom(200, 1, 0.5)  # Mengganti x2 dengan status_langganan (0 = Baru, 1 = Setia)
frekuensi_kunjungan <- rnorm(200)     # Mengganti x3 dengan frekuensi_kunjungan

# Menghitung linier prediktor berdasarkan koefisien yang sesuai dengan kasus
lin_pred <- -0.5 + 1.2 * usia - 0.8 * status_langganan + 0.5 * frekuensi_kunjungan  # Koefisien disesuaikan dari analisis sebelumnya
p <- 1 / (1 + exp(-lin_pred))
pembelian <- rbinom(200, 1, p)        # Mengganti y dengan pembelian

# Membuat data frame sesuai kasus
df_pemasaran <- data.frame(pembelian = pembelian, usia, status_langganan, frekuensi_kunjungan)

# Membangun model regresi logistik
model <- glm(pembelian ~ usia + status_langganan + frekuensi_kunjungan, data = df_pemasaran, family = binomial)

# Mendapatkan probabilitas prediksi
prob <- predict(model, type = "response")

# Membuat kurva Precision-Recall
# scores.class0: probabilitas untuk kelas positif (pembelian = 1)
# scores.class1: probabilitas untuk kelas negatif (pembelian = 0)
pr <- pr.curve(scores.class0 = prob[df_pemasaran$pembelian == 1],
               scores.class1 = prob[df_pemasaran$pembelian == 0],
               curve = TRUE)

# Plot kurva Precision-Recall
plot(pr)

dimana :

  • Sumbu X (Recall) : Recall bervariasi dari 0.0 (tidak ada pembeli terdeteksi) hingga 1.0 (semua pembeli terdeteksi), sesuai dengan perubahan threshold dari tinggi ke rendah

  • Sumbu Y (Precision) : Precision mendekati 0.8 pada recall rendah, menurun seiring recall meningkat, yang menunjukkan tradeoff tipikal.

  • Area Under Cover : Nilai 0.6819579 menunjukkan performa yang cukup baik

  • Warna pada Kurva : Gradien warna (dari ungu hingga merah) kemungkinan menunjukkan perubahan threshold atau distribusi data di sepanjang kurva, dengan ungu mewakili threshold tinggi (precision tinggi, recall rendah) dan merah mewakili threshold rendah (recall tinggi, precision rendah).

Kesimpulan :

  • Kurva dimulai di titik tinggi pada sumbu Y (precision ~0.8) dengan recall rendah (~0.0), yang menunjukkan bahwa pada threshold tinggi (misalnya, 0.90), model sangat selektif dan akurat dalam memprediksi pembeli, tetapi hanya mendeteksi sedikit dari mereka. Seiring recall meningkat (threshold menurun, misalnya dari 0.90 ke 0.10), precision menurun secara bertahap hingga mendekati 0.2 pada recall ~1.0. Ini mencerminkan peningkatan false positives saat model berusaha mendeteksi lebih banyak pembeli.

  • Recall Tinggi (0.8 - 1.0): Precision turun di bawah 0.3, menunjukkan bahwa pada threshold rendah (misalnya, 0.15 dengan sensitivitas 0.9275), banyak false positives muncul, yang kurang efisien untuk targeting.

  • Recall Sedang (0.4 - 0.6): Precision tetap di atas 0.4 - 0.5, yang sesuai dengan threshold 0.35–0.50 (sensitivitas 0.6957–0.5362), memberikan keseimbangan yang baik berdasarkan matriks kebingungan sebelumnya (precision 0.7115).
  • Recall Rendah (0.0 - 0.2): Precision tinggi (~0.7 - 0.8), sesuai dengan threshold tinggi (misalnya, 0.70 dengan sensitivitas 0.3043), yang ideal untuk precision maksimum.

Distribusi Multinomial

Distribusi multinomial menggambarkan distribusi probabilitas dari hasil \(n\) percobaan independen, di mana setiap percobaan dapat menghasilkan salah satu dari \(k\) kemungkinan hasil (kategori), dengan probabilitas masing-masing kategori tetap. Misalnya, jika Anda melempar dadu enam sisi \(n\) kali, distribusi multinomial dapat digunakan untuk memodelkan jumlah kemunculan setiap angka (1, 2, 3, 4, 5, 6).

Parameter

Fungsi Probabilitas

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

Di mana :

Penjelasan :

Fraksi \(\frac{n!}{x_1! \cdot x_2! \cdot \ldots \cdot x_k!}\) adalah koefisien multinomial, yang menghitung jumlah cara untuk membagi \(n\) percobaan ke dalam \(k\) kategori dengan frekuensi \(x_1, x_2, \ldots, x_k\).

\(p_1^{x_1} \cdot p_2^{x_2} \cdot \ldots \cdot p_k^{x_k}\) adalah probabilitas bersama dari hasil tersebut, yang diasumsikan independen di setiap percobaan.

1. Sifat Distribusi Multinomial

Distribusi multinomial memiliki beberapa sifat penting yang relevan untuk analisis data kategorikal:


2. Aplikasi Distribusi Multinomial dalam Buku

Distribusi multinomial adalah dasar untuk banyak metode analisis data kategorikal:


3. Inferensi untuk Distribusi Multinomial

Inferensi statistik untuk parameter distribusi multinomial (\(\pi_1, \pi_2, \ldots, \pi_k\)) melibatkan estimasi dan pengujian hipotesis:


4. Hubungan dengan Distribusi Lain

Contoh Kasus

Dalam era digital, perusahaan teknologi seperti xAI sering kali mengandalkan analisis sentimen media sosial untuk memahami persepsi pengguna terhadap produk atau layanan mereka, misalnya, asisten AI seperti Grok. Studi kasus ini berfokus pada analisis sentimen dari postingan di platform X terkait fitur baru Grok 3 yang diluncurkan pada Juni 2025. Tujuannya adalah untuk mengklasifikasikan sentimen pengguna ke dalam tiga kategori: positif, negatif, dan netral, berdasarkan teks postingan, dan menggunakan distribusi multinomial untuk memodelkan distribusi sentimen.

Hasil analisis sentimen :

  • Positif : 350 pengguna

  • Netral : 400 pengguna

  • Negatif : 250 pengguna

Probabilitas :

  • \(P_p\) : 0.35

  • \(P_Net\) : 0.4

  • \(P_Neg\) : 0.25

Ditanya : Berapa peluang pada 1000 pengguna akan ada 350 pengguna berkomentar postif, 400 pengguna berkomentar netral dan 250 orang berkomentar negatif?

n <- 1000
x <- c(350, 400, 250)

# Hitung proporsi aktual dari data
p <- x / n  # [0.35, 0.40, 0.25]

# Hitung koefisien menggunakan log untuk menghindari overflow
log_faktorial_total <- lfactorial(n)
log_faktorial_x <- sum(lfactorial(x))
log_koefisien <- log_faktorial_total - log_faktorial_x
koefisien <- exp(log_koefisien)

# Hitung peluang
log_peluang <- log_koefisien + sum(x * log(p))
peluang <- exp(log_peluang)

print(peluang)
## [1] 0.0008501267

Interpretasi :

Probabilitas bahwa 1000 pengguna akan memberikan 350 komentar posisitif, 400 komentar netral, dan 250 komentar negatif adalah 0.0008501267

Multinomial Regresi Logistik

Regresi logistik multinomial digunakan ketika variabel dependen memiliki \(k\) kategori (di mana \(k \geq 3\)) dan tidak ada urutan alami di antara kategori-kategori tersebut. Berbeda dengan regresi logistik biner yang memodelkan peluang log (log-odds) untuk dua hasil, regresi multinomial memodelkan hubungan antara variabel independen (prediktor) dan peluang relatif untuk setiap kategori dibandingkan dengan kategori referensi (baseline).

Baseline Category Logit Model

Model baseline category logit adalah pendekatan utama dalam regresi logistik multinomial. Model ini membandingkan peluang log (log-odds) setiap kategori terhadap satu kategori referensi (baseline), yang dipilih sebagai dasar perbandingan.

Misalkan \(Y\) adalah variabel dependen dengan \(k\) kategori, dan kategori ke-\(j\) dipilih sebagai baseline. Untuk setiap kategori lain \(i\) (di mana \(i \neq j\)), model menghitung log-odds rasio peluang kategori \(i\) terhadap kategori \(j\).

Rumus Model: Untuk kategori \(i\) (dibandingkan dengan baseline \(j\)):

\[ \log\left(\frac{P(Y = i | X)}{P(Y = j | X)}\right) = \beta_{0i} + \beta_{1i}X_1 + \beta_{2i}X_2 + \ldots + \beta_{pi}X_p \]

Di mana: - \(P(Y = i | X)\): Probabilitas kategori \(i\) diberikan prediktor \(X\), - \(P(Y = j | X)\): Probabilitas baseline kategori \(j\), - \(\beta_{0i}\): Intersep untuk kategori \(i\) vs. \(j\), - \(\beta_{1i}, \beta_{2i}, \ldots, \beta_{pi}\): Koefisien untuk prediktor \(X_1, X_2, \ldots, X_p\) pada perbandingan \(i\) vs. \(j\).

Probabilitas: Probabilitas untuk setiap kategori dihitung menggunakan softmax:

\[ P(Y = i | X) = \frac{\exp(\beta_{0i} + \beta_{1i}X_1 + \ldots + \beta_{pi}X_p)}{\sum_{m=1}^{k} \exp(\beta_{0m} + \beta_{1m}X_1 + \ldots + \beta_{pm}X_p)} \]

Di mana penyebut adalah normalisasi untuk memastikan jumlah probabilitas semua kategori adalah 1.

Estimasi Parameter

Parameter dalam regresi logistik multinomial diestimasi menggunakan metode maksimum likelihood estimation (MLE), yang memaksimalkan fungsi likelihood berdasarkan data observasi. Prosesnya melibatkan:

  1. Fungsi Likelihood:

    • Untuk setiap observasi, likelihood adalah probabilitas kategori yang diamati berdasarkan model.

    • Fungsi likelihood keseluruhan adalah produk likelihood untuk semua observasi.

Contoh Kasus

Sebuah perusahaan logistik ingin memahami faktor-faktor yang memengaruhi preferensi karyawan terhadap jenis alat transportasi kerja yang disediakan: mobil, sepeda motor, atau sepeda. Tujuannya adalah untuk mengoptimalkan kebijakan transportasi guna meningkatkan efisiensi dan kepuasan karyawan di berbagai rute pengiriman pada Juni 2025.

  1. Data : Perusahaan melakukan survei terhadap 150 karyawan kurir dan mengumpulkan data berikut:
    • Transport: Jenis alat transportasi yang dipilih (mobil, sepeda motor, sepeda).

    • Age: Usia karyawan.

    • Department : Departemen tempat kurir tersebut bekerja (Logistik, Operation, Admin)

    • Experience: Tahun pengalaman kerja.

Input Data

set.seed(36)
n <- 150
Department <- sample(c("Logistics", "Operations", "Admin"), n, replace = TRUE)  # Mengganti dengan departemen logistik
Age <- round(rnorm(n, mean = 35, sd = 7))  # Usia karyawan
Experience <- round(pmax(rnorm(n, mean = 7, sd = 3), 0))  # Pengalaman kerja (minimal 0)

# Simulasikan Transport berdasarkan probabilitas berbeda per Department
Transport <- sapply(Department, function(dep) {
  if (dep == "Logistics") {
    sample(c("Mobil", "Sepeda Motor", "Sepeda"), size = 1, prob = c(0.6, 0.3, 0.1))  # Preferensi untuk rute berat
  } else if (dep == "Operations") {
    sample(c("Mobil", "Sepeda Motor", "Sepeda"), size = 1, prob = c(0.3, 0.5, 0.2))  # Preferensi untuk rute sedang
  } else {
    sample(c("Mobil", "Sepeda Motor", "Sepeda"), size = 1, prob = c(0.2, 0.3, 0.5))  # Preferensi untuk rute ringan
  }
})

df <- data.frame(Transport = factor(Transport), Age, Department = factor(Department), Experience)
df$Transport <- relevel(df$Transport, ref = "Sepeda")  # Baseline: Sepeda

# Tampilkan 6 baris pertama
head(df)

Estimasi Model

library(nnet)
model_mnlogit <- multinom(Transport ~ Age + Department + Experience, data = df)
## # weights:  18 (10 variable)
## initial  value 164.791843 
## iter  10 value 138.498151
## final  value 137.898479 
## converged
# Ringkasan hasil
summary(model_mnlogit)
## Call:
## multinom(formula = Transport ~ Age + Department + Experience, 
##     data = df)
## 
## Coefficients:
##              (Intercept)         Age DepartmentLogistics DepartmentOperations
## Mobil          -5.282244 0.103026832            3.544800             1.768236
## Sepeda Motor   -1.512575 0.009075642            2.446303             1.634146
##              Experience
## Mobil        0.05121187
## Sepeda Motor 0.05426996
## 
## Std. Errors:
##              (Intercept)        Age DepartmentLogistics DepartmentOperations
## Mobil           1.637725 0.03880227           0.7493704            0.5761138
## Sepeda Motor    1.397679 0.03419644           0.7304759            0.5063352
##              Experience
## Mobil        0.08448080
## Sepeda Motor 0.07850996
## 
## Residual Deviance: 275.797 
## AIC: 295.797

Interpretasi :

  • Departemen: Departemen Logistics dan Operations secara signifikan memengaruhi preferensi mobil dan sepeda motor dibandingkan sepeda, dengan Logistics menunjukkan pengaruh terkuat (odds ratio tinggi).
  • Usia: Hanya signifikan untuk mobil, menunjukkan karyawan yang lebih tua cenderung memilih mobil, mungkin karena kenyamanan.
  • Pengalaman: Tidak signifikan, menunjukkan pengalaman kerja tidak memengaruhi pilihan transportasi secara nyata.
z <- summary(model_mnlogit)$coefficients / summary(model_mnlogit)$standard.errors
pval <- 2 * (1 - pnorm(abs(z)))
round(pval, 4)
##              (Intercept)    Age DepartmentLogistics DepartmentOperations
## Mobil             0.0013 0.0079               0e+00               0.0021
## Sepeda Motor      0.2792 0.7907               8e-04               0.0012
##              Experience
## Mobil            0.5444
## Sepeda Motor     0.4894

Kesimpulan :

  • Usia: Hanya signifikan untuk “Mobil” vs. “Sepeda”, menunjukkan preferensi mobil meningkat seiring usia.
  • Departemen: Baik “Logistics” maupun “Operations” signifikan untuk kedua kategori (“Mobil” dan “Sepeda Motor”) vs. “Sepeda”, dengan efek terkuat dari Logistics pada mobil.
  • Intercept: Signifikan hanya untuk “Mobil”, mencerminkan perbedaan peluang dasar yang besar.

  • Pengalaman: Tidak ada efek nyata pada pilihan transportasi, mungkin karena preferensi lebih dipengaruhi oleh kebutuhan kerja daripada lama kerja.

  • Intercept untuk Sepeda Motor: Tidak signifikan, menunjukkan bahwa perbedaan peluang dasar dengan sepeda kecil tanpa prediktor.

  • Logistics: Karyawan di departemen ini cenderung memilih mobil (kemungkinan karena rute berat) dan sepeda motor, sesuai dengan probabilitas simulasi (60% mobil, 30% sepeda motor).

  • Operations: Lebih condong ke sepeda motor (50% dalam simulasi), yang konsisten dengan efek signifikan.
  • Admin: Sebagai baseline, lebih banyak memilih sepeda (50% dalam simulasi), yang mendukung model.

Prediksi dan Validasi

df$Predicted <- predict(model_mnlogit)
table(Predicted = df$Predicted, Actual = df$Transport)
##               Actual
## Predicted      Sepeda Mobil Sepeda Motor
##   Sepeda           24     7           11
##   Mobil             4    37           20
##   Sepeda Motor     12    13           22

Kesimpulan Akhir : Model memiliki akurasi sedang (sekitar 55%), menunjukkan performa yang tidak optimal namun masih berguna untuk analisis awal.

Visualisasi

library(ggplot2)

# Asumsi df sudah berisi kolom Transport, Age, Department, Experience, dan Predicted
ggplot(df, aes(x = Age, y = Experience, color = Predicted)) +
  geom_point(size = 2) +
  labs(title = "Multinomial Logistic Regression Predictions",
       x = "Age", y = "Experience") +
  theme_minimal()

Kesimpulan Visualisasi :

Berdasarkan visualisasi plot scatter dari model regresi logistik multinomial yang menunjukkan prediksi preferensi alat transportasi kerja (Sepeda, Mobil, Sepeda Motor) berdasarkan usia dan pengalaman karyawan, terlihat bahwa distribusi prediksi cukup bervariasi di seluruh rentang usia (20-50 tahun) dan pengalaman (0-10 tahun). Karyawan yang lebih tua (usia > 40) cenderung diprediksi memilih Mobil (hijau), yang konsisten dengan efek signifikan usia pada prediksi Mobil dari analisis sebelumnya, sementara Sepeda (merah) lebih dominan pada usia muda (< 30) dan pengalaman rendah, sesuai dengan baseline. Sepeda Motor (biru) tersebar merata, namun tampaknya lebih sering diprediksi untuk usia dan pengalaman menengah, mencerminkan performa model yang lebih lemah untuk kategori ini seperti yang terlihat pada matriks kebingungan. Secara keseluruhan, plot ini mendukung strategi alokasi transportasi berdasarkan usia dan departemen, meskipun perlu penyempurnaan model untuk meningkatkan akurasi prediksi, terutama untuk Sepeda Motor.

Regresi Logistik Ordinal

Regresi logistik ordinal dirancang untuk memprediksi hubungan antara variabel independen (prediktor) dan variabel dependen ordinal. Model ini mengasumsikan bahwa hubungan antara prediktor dan log-odds kategori bertambah secara konsisten sepanjang skala ordinal. Hal ini membuatnya lebih efisien dibandingkan regresi multinomial, terutama ketika jumlah kategori meningkat.

Kumulatif Logit Model

Model kumulatif logit adalah pendekatan utama dalam regresi logistik ordinal. Model ini memodelkan log-odds kumulatif dari kategori yang lebih rendah atau sama terhadap kategori yang lebih tinggi, menggunakan ambang batas (cutpoints) untuk memisahkan kategori.

Definisi: Misalkan \(Y\) adalah variabel dependen ordinal dengan \(k\) kategori (misalnya, 1, 2, 3). Definisikan \(\pi_j = P(Y \leq j | X)\) sebagai probabilitas kumulatif bahwa \(Y\) berada pada kategori \(j\) atau lebih rendah, diberikan prediktor \(X\).

Log-odds kumulatif untuk kategori \(j\) vs. kategori di atasnya adalah:

\[ \log\left(\frac{P(Y \leq j | X)}{P(Y > j | X)}\right) = \log\left(\frac{\pi_j}{1 - \pi_j}\right) \]

Rumus Model: Model kumulatif logit menyatakan bahwa log-odds ini adalah fungsi linier dari prediktor:

\[ \log\left(\frac{P(Y \leq j | X)}{P(Y > j | X)}\right) = \alpha_j - (\beta_1 X_1 + \beta_2 X_2 + \ldots + \beta_p X_p) \]

Di mana:

- \(\alpha_j\): Ambang batas (cutpoint) untuk kategori \(j\) (berbeda untuk setiap \(j\), tetapi jumlahnya \(k-1\)).

- \(\beta_1, \beta_2, \ldots, \beta_p\): Koefisien untuk prediktor \(X_1, X_2, \ldots, X_p\), yang sama di seluruh ambang batas (proportional odds).

- \(P(Y > j | X) = 1 - \pi_j\).

Probabilitas Kategori: Probabilitas untuk setiap kategori \(j\) dihitung sebagai:

\[ P(Y = j | X) = P(Y \leq j | X) - P(Y \leq j-1 | X) \]

Di mana \(P(Y \leq 0 | X) = 0\) dan \(P(Y \leq k | X) = 1\).

1. Asumsi Rasio Odds Proporsional

Model ini mengasumsikan bahwa rasio odds untuk membandingkan \(Y \leq j\) dengan \(Y > j\) adalah konstan untuk semua \(j\), hanya bergeser oleh \(\alpha_j\). Rasio odds untuk perubahan satu unit pada prediktor \(x_k\) adalah: \[ \text{Odds Ratio} = \exp(\beta_k) \] Ini berarti efek prediktor konsisten di semua batas kategori.

Alternatif Fungsi Link

Selain logit kumulatif, model lain seperti probit kumulatif (berdasarkan distribusi normal) atau log-log kumulatif dapat digunakan, tetapi logit lebih umum karena interpretasi rasio odds yang intuitif.

Contoh

Dalam data kepuasan kerja, model logit kumulatif menunjukkan bahwa kenaikan pendapatan meningkatkan odds untuk berada pada kategori kepuasan yang lebih tinggi, dengan koefisien \(\beta\) yang diinterpretasikan sebagai log rasio odds.


2. Interpretasi Parameter

Interpretasi parameter dalam regresi logistik ordinal berfokus pada efek prediktor terhadap probabilitas kumulatif:

  • Koefisien \(\beta_k\):
    Nilai \(\beta_k > 0\) menunjukkan bahwa kenaikan pada prediktor \(x_k\) meningkatkan probabilitas berada pada kategori respons yang lebih tinggi (atau mengurangi probabilitas \(Y \leq j\)). Sebaliknya, \(\beta_k < 0\) menunjukkan efek sebaliknya.
    Rasio odds \(\exp(\beta_k)\) mengukur perubahan odds untuk \(Y \leq j\) vs. \(Y > j\) per unit kenaikan \(x_k\).

  • Intersep \(\alpha_j\):
    Parameter \(\alpha_j\) menentukan batas antara kategori \(j\) dan \(j+1\) ketika semua prediktor \(\mathbf{x} = 0\). Nilai ini tidak selalu memiliki interpretasi praktis kecuali prediktor distandarisasi.

  • Probabilitas Kategori:
    Probabilitas untuk kategori spesifik \(j\) dapat dihitung sebagai: \[ P(Y = j | \mathbf{x}) = P(Y \leq j | \mathbf{x}) - P(Y \leq j-1 | \mathbf{x}) \] dengan: \[ P(Y \leq j | \mathbf{x}) = \frac{\exp(\alpha_j - \boldsymbol{\beta}^T \mathbf{x})}{1 + \exp(\alpha_j - \boldsymbol{\beta}^T \mathbf{x})} \]

Contoh

Dalam studi kepuasan kerja, koefisien \(\beta = 0.5\) untuk pendapatan menunjukkan bahwa setiap kenaikan \(1000\) dalam pendapatan meningkatkan odds untuk kepuasan yang lebih tinggi sebesar \(\exp(0.5) \approx 1.65\).


3. Inferensi untuk Regresi Logistik Ordinal

Inferensi statistik melibatkan estimasi parameter, pengujian signifikansi, dan pemeriksaan asumsi model:

  • Estimasi Parameter:
    Parameter \(\alpha_j\) dan \(\boldsymbol{\beta}\) diestimasi menggunakan Maximum Likelihood Estimation (MLE), dengan distribusi multinomial untuk data. Perangkat lunak seperti SAS atau R digunakan untuk perhitungan.

  • Interval Kepercayaan:
    Interval kepercayaan untuk \(\beta_k\) dihitung menggunakan pendekatan Wald: \[ \hat{\beta}_k \pm z_{\alpha/2} \cdot SE(\hat{\beta}_k) \] Untuk rasio odds, interval kepercayaan dihitung dengan mengambil eksponen dari batas-batas tersebut.

  • Pengujian Signifikansi:

    • Uji Wald: Menguji hipotesis nol \(H_0: \beta_k = 0\) dengan statistik: \[ z = \frac{\hat{\beta}_k}{SE(\hat{\beta}_k)} \]
    • Uji Rasio Likelihud: Membandingkan model dengan dan tanpa prediktor tertentu untuk mengevaluasi signifikansi.
  • Goodness-of-Fit:
    Kecocokan model dapat dievaluasi menggunakan statistik deviansi (\(G^2\)) atau uji seperti Hosmer-Lemeshow untuk data ordinal. Namun, uji ini lebih kompleks untuk model ordinal dibandingkan model biner.

  • Pemeriksaan Asumsi Rasio Odds Proporsional:
    Asumsi bahwa \(\boldsymbol{\beta}\) konstan untuk semua \(j\) dapat diuji menggunakan score test atau Brant test. Jika asumsi ini dilanggar, model alternatif seperti partial proportional odds model atau model multinomial dapat digunakan.


4. Model Alternatif untuk Data Ordinal

Jika asumsi rasio odds proporsional tidak terpenuhi, beberapa model alternatif dapat dipertimbangkan:

  • Model Logit Berdekatan (Adjacent-Categories Logit):
    Memodelkan log rasio odds untuk kategori berdekatan: \[ \log\left(\frac{P(Y = j | \mathbf{x})}{P(Y = j+1 | \mathbf{x})}\right) = \alpha_j - \boldsymbol{\beta}^T \mathbf{x} \] Model ini cocok untuk data dengan fokus pada perbandingan kategori berdekatan.

  • Model Logit Kontinuasi (Continuation-Ratio Logit):
    Memodelkan odds untuk mencapai kategori tertentu dibandingkan tetap di kategori sebelumnya: \[ \log\left(\frac{P(Y = j | \mathbf{x})}{P(Y < j | \mathbf{x})}\right) = \alpha_j - \boldsymbol{\beta}^T \mathbf{x} \] Model ini relevan untuk proses berurutan, seperti tingkat pendidikan.

  • Model Probit atau Log-Log:
    Menggunakan fungsi link probit atau log-log untuk menangani distribusi yang berbeda, meskipun interpretasinya kurang intuitif dibandingkan logit.

Contoh Kasus

set.seed(36)
n <- 200
Waiting_Time <- round(runif(n, 1, 10))  # Mengganti speed dengan Waiting_Time (durasi tunggu dalam menit)
Satisfaction <- cut(5 + 0.5 * Waiting_Time + rnorm(n),  # Mengganti satisfaction dengan Satisfaction
                   breaks = c(-Inf, 5.5, 7.5, Inf),    # Ambang batas untuk kategori ordinal
                   labels = c("Tidak Puas", "Cukup", "Puas"),  # Label kategori
                   ordered_result = TRUE)
df <- data.frame(Satisfaction, Waiting_Time)  # Mengganti df dengan kolom yang sesuai
head(df)

Estimasi Model Ordinal

library(MASS)
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
model_ord <- polr(Satisfaction ~ Waiting_Time, data = df, Hess = TRUE)
summary(model_ord)
## Call:
## polr(formula = Satisfaction ~ Waiting_Time, data = df, Hess = TRUE)
## 
## Coefficients:
##               Value Std. Error t value
## Waiting_Time 0.9068     0.1095   8.278
## 
## Intercepts:
##                  Value  Std. Error t value
## Tidak Puas|Cukup 0.5994 0.4124     1.4534 
## Cukup|Puas       4.2008 0.5707     7.3612 
## 
## Residual Deviance: 208.9792 
## AIC: 214.9792

Kesimpulan :

  • Koefisien 0.9068 menunjukkan bahwa untuk setiap peningkatan satu menit dalam durasi tunggu, log-odds (peluang relatif) pasien berada di kategori kepuasan yang lebih tinggi (misalnya, dari “Tidak Puas” ke “Cukup” atau dari “Cukup” ke “Puas”) meningkat sebesar 0.9068, dengan asumsi variabel lain konstan.

  • Durasi tunggu (Waiting_Time) memiliki efek positif yang signifikan terhadap kepuasan pasien. Artinya, semakin lama pasien menunggu (dalam rentang 1-10 menit), semakin tinggi peluang mereka untuk berada di kategori kepuasan yang lebih tinggi (misalnya, “Puas” dibandingkan “Cukup” atau “Tidak Puas”).

Nilai P-Value

(ctable <- coef(summary(model_ord)))
##                      Value Std. Error  t value
## Waiting_Time     0.9068327  0.1095421 8.278394
## Tidak Puas|Cukup 0.5993938  0.4124128 1.453383
## Cukup|Puas       4.2007859  0.5706651 7.361210
p <- pnorm(abs(ctable[, "t value"]), lower.tail = FALSE) * 2
(ctable <- cbind(ctable, "p value" = round(p, 4)))
##                      Value Std. Error  t value p value
## Waiting_Time     0.9068327  0.1095421 8.278394  0.0000
## Tidak Puas|Cukup 0.5993938  0.4124128 1.453383  0.1461
## Cukup|Puas       4.2007859  0.5706651 7.361210  0.0000

Kesimpulan :

  • Nilai p-value sebesar 0.0000 (sangat kecil, jauh di bawah ambang batas 0.05) menunjukkan bahwa variabel Waiting_Time (durasi tunggu) memiliki efek yang sangat signifikan secara statistik terhadap tingkat kepuasan pasien. Sehingga Durasi tunggu memiliki pengaruh yang sangat signifikan terhadap kepuasan pasien. Rumah sakit dapat mempercayai bahwa durasi tunggu adalah faktor penting dalam menentukan tingkat kepuasan pasien di UGD

  • Nilai p-value sebesar 0.1461 (lebih besar dari 0.05) menunjukkan bahwa ambang batas (intersep) antara kategori “Tidak Puas” dan “Cukup” tidak signifikan secara statistik pada tingkat kepercayaan 95% ($\alpha$ = 0.05).

  • Nilai p-value sebesar 0.0000 (sangat kecil, jauh di bawah 0.05) menunjukkan bahwa ambang batas antara kategori “Cukup” dan “Puas” sangat signifikan secara statistik. Sehingga hal ini menunjukkan bahwa faktor durasi tunggu lebih kuat dalam memengaruhi transisi ke kategori kepuasan yang lebih tinggi

Prediksi Probabilitas

newdata <- data.frame(speed = 5:9)
predict(model_ord, newdata = newdata, type = "probs")
## Warning: 'newdata' had 5 rows but variables found have 200 rows
##       Tidak Puas       Cukup       Puas
## 1   0.0031776352 0.101430219 0.89539215
## 2   0.0031776352 0.101430219 0.89539215
## 3   0.0012855683 0.043764868 0.95494956
## 4   0.1070624427 0.707553600 0.18538396
## 5   0.0012855683 0.043764868 0.95494956
## 6   0.4237399822 0.540480805 0.03577921
## 7   0.0002098446 0.007423691 0.99236646
## 8   0.0191753167 0.398243691 0.58258099
## 9   0.2289460980 0.686893791 0.08416011
## 10  0.0461796049 0.593380420 0.36043997
## 11  0.0461796049 0.593380420 0.36043997
## 12  0.0078325591 0.216566804 0.77560064
## 13  0.0078325591 0.216566804 0.77560064
## 14  0.4237399822 0.540480805 0.03577921
## 15  0.2289460980 0.686893791 0.08416011
## 16  0.0191753167 0.398243691 0.58258099
## 17  0.0078325591 0.216566804 0.77560064
## 18  0.2289460980 0.686893791 0.08416011
## 19  0.0078325591 0.216566804 0.77560064
## 20  0.1070624427 0.707553600 0.18538396
## 21  0.0461796049 0.593380420 0.36043997
## 22  0.0461796049 0.593380420 0.36043997
## 23  0.0461796049 0.593380420 0.36043997
## 24  0.0012855683 0.043764868 0.95494956
## 25  0.0031776352 0.101430219 0.89539215
## 26  0.0461796049 0.593380420 0.36043997
## 27  0.0005195121 0.018173995 0.98130649
## 28  0.0012855683 0.043764868 0.95494956
## 29  0.0078325591 0.216566804 0.77560064
## 30  0.0191753167 0.398243691 0.58258099
## 31  0.0078325591 0.216566804 0.77560064
## 32  0.0031776352 0.101430219 0.89539215
## 33  0.0005195121 0.018173995 0.98130649
## 34  0.1070624427 0.707553600 0.18538396
## 35  0.0461796049 0.593380420 0.36043997
## 36  0.0078325591 0.216566804 0.77560064
## 37  0.0031776352 0.101430219 0.89539215
## 38  0.1070624427 0.707553600 0.18538396
## 39  0.0461796049 0.593380420 0.36043997
## 40  0.0002098446 0.007423691 0.99236646
## 41  0.0005195121 0.018173995 0.98130649
## 42  0.0012855683 0.043764868 0.95494956
## 43  0.2289460980 0.686893791 0.08416011
## 44  0.0078325591 0.216566804 0.77560064
## 45  0.0005195121 0.018173995 0.98130649
## 46  0.0031776352 0.101430219 0.89539215
## 47  0.2289460980 0.686893791 0.08416011
## 48  0.0461796049 0.593380420 0.36043997
## 49  0.1070624427 0.707553600 0.18538396
## 50  0.0078325591 0.216566804 0.77560064
## 51  0.0002098446 0.007423691 0.99236646
## 52  0.0012855683 0.043764868 0.95494956
## 53  0.4237399822 0.540480805 0.03577921
## 54  0.0002098446 0.007423691 0.99236646
## 55  0.0005195121 0.018173995 0.98130649
## 56  0.2289460980 0.686893791 0.08416011
## 57  0.0012855683 0.043764868 0.95494956
## 58  0.0002098446 0.007423691 0.99236646
## 59  0.0461796049 0.593380420 0.36043997
## 60  0.0012855683 0.043764868 0.95494956
## 61  0.0461796049 0.593380420 0.36043997
## 62  0.0031776352 0.101430219 0.89539215
## 63  0.1070624427 0.707553600 0.18538396
## 64  0.2289460980 0.686893791 0.08416011
## 65  0.0005195121 0.018173995 0.98130649
## 66  0.2289460980 0.686893791 0.08416011
## 67  0.1070624427 0.707553600 0.18538396
## 68  0.0461796049 0.593380420 0.36043997
## 69  0.0012855683 0.043764868 0.95494956
## 70  0.2289460980 0.686893791 0.08416011
## 71  0.2289460980 0.686893791 0.08416011
## 72  0.2289460980 0.686893791 0.08416011
## 73  0.0005195121 0.018173995 0.98130649
## 74  0.0005195121 0.018173995 0.98130649
## 75  0.0078325591 0.216566804 0.77560064
## 76  0.0191753167 0.398243691 0.58258099
## 77  0.0012855683 0.043764868 0.95494956
## 78  0.4237399822 0.540480805 0.03577921
## 79  0.2289460980 0.686893791 0.08416011
## 80  0.0012855683 0.043764868 0.95494956
## 81  0.0012855683 0.043764868 0.95494956
## 82  0.0031776352 0.101430219 0.89539215
## 83  0.0461796049 0.593380420 0.36043997
## 84  0.0012855683 0.043764868 0.95494956
## 85  0.0461796049 0.593380420 0.36043997
## 86  0.0078325591 0.216566804 0.77560064
## 87  0.0012855683 0.043764868 0.95494956
## 88  0.4237399822 0.540480805 0.03577921
## 89  0.0078325591 0.216566804 0.77560064
## 90  0.0191753167 0.398243691 0.58258099
## 91  0.0031776352 0.101430219 0.89539215
## 92  0.0012855683 0.043764868 0.95494956
## 93  0.0078325591 0.216566804 0.77560064
## 94  0.2289460980 0.686893791 0.08416011
## 95  0.0002098446 0.007423691 0.99236646
## 96  0.0005195121 0.018173995 0.98130649
## 97  0.0031776352 0.101430219 0.89539215
## 98  0.0461796049 0.593380420 0.36043997
## 99  0.0078325591 0.216566804 0.77560064
## 100 0.0031776352 0.101430219 0.89539215
## 101 0.0005195121 0.018173995 0.98130649
## 102 0.0002098446 0.007423691 0.99236646
## 103 0.0078325591 0.216566804 0.77560064
## 104 0.0078325591 0.216566804 0.77560064
## 105 0.0078325591 0.216566804 0.77560064
## 106 0.0078325591 0.216566804 0.77560064
## 107 0.0031776352 0.101430219 0.89539215
## 108 0.0031776352 0.101430219 0.89539215
## 109 0.1070624427 0.707553600 0.18538396
## 110 0.1070624427 0.707553600 0.18538396
## 111 0.0078325591 0.216566804 0.77560064
## 112 0.0002098446 0.007423691 0.99236646
## 113 0.0031776352 0.101430219 0.89539215
## 114 0.0078325591 0.216566804 0.77560064
## 115 0.4237399822 0.540480805 0.03577921
## 116 0.4237399822 0.540480805 0.03577921
## 117 0.0461796049 0.593380420 0.36043997
## 118 0.0191753167 0.398243691 0.58258099
## 119 0.1070624427 0.707553600 0.18538396
## 120 0.0078325591 0.216566804 0.77560064
## 121 0.0031776352 0.101430219 0.89539215
## 122 0.0191753167 0.398243691 0.58258099
## 123 0.0012855683 0.043764868 0.95494956
## 124 0.0078325591 0.216566804 0.77560064
## 125 0.0461796049 0.593380420 0.36043997
## 126 0.1070624427 0.707553600 0.18538396
## 127 0.0078325591 0.216566804 0.77560064
## 128 0.0012855683 0.043764868 0.95494956
## 129 0.2289460980 0.686893791 0.08416011
## 130 0.0191753167 0.398243691 0.58258099
## 131 0.0012855683 0.043764868 0.95494956
## 132 0.2289460980 0.686893791 0.08416011
## 133 0.0461796049 0.593380420 0.36043997
## 134 0.0078325591 0.216566804 0.77560064
## 135 0.0012855683 0.043764868 0.95494956
## 136 0.0005195121 0.018173995 0.98130649
## 137 0.0012855683 0.043764868 0.95494956
## 138 0.0002098446 0.007423691 0.99236646
## 139 0.1070624427 0.707553600 0.18538396
## 140 0.0002098446 0.007423691 0.99236646
## 141 0.0031776352 0.101430219 0.89539215
## 142 0.2289460980 0.686893791 0.08416011
## 143 0.0191753167 0.398243691 0.58258099
## 144 0.0461796049 0.593380420 0.36043997
## 145 0.0191753167 0.398243691 0.58258099
## 146 0.2289460980 0.686893791 0.08416011
## 147 0.1070624427 0.707553600 0.18538396
## 148 0.0005195121 0.018173995 0.98130649
## 149 0.0461796049 0.593380420 0.36043997
## 150 0.1070624427 0.707553600 0.18538396
## 151 0.1070624427 0.707553600 0.18538396
## 152 0.0031776352 0.101430219 0.89539215
## 153 0.0078325591 0.216566804 0.77560064
## 154 0.0005195121 0.018173995 0.98130649
## 155 0.0078325591 0.216566804 0.77560064
## 156 0.0078325591 0.216566804 0.77560064
## 157 0.0078325591 0.216566804 0.77560064
## 158 0.4237399822 0.540480805 0.03577921
## 159 0.0191753167 0.398243691 0.58258099
## 160 0.0078325591 0.216566804 0.77560064
## 161 0.2289460980 0.686893791 0.08416011
## 162 0.0012855683 0.043764868 0.95494956
## 163 0.0005195121 0.018173995 0.98130649
## 164 0.0191753167 0.398243691 0.58258099
## 165 0.0012855683 0.043764868 0.95494956
## 166 0.4237399822 0.540480805 0.03577921
## 167 0.0012855683 0.043764868 0.95494956
## 168 0.0031776352 0.101430219 0.89539215
## 169 0.0005195121 0.018173995 0.98130649
## 170 0.0002098446 0.007423691 0.99236646
## 171 0.0078325591 0.216566804 0.77560064
## 172 0.0078325591 0.216566804 0.77560064
## 173 0.0012855683 0.043764868 0.95494956
## 174 0.0078325591 0.216566804 0.77560064
## 175 0.4237399822 0.540480805 0.03577921
## 176 0.0191753167 0.398243691 0.58258099
## 177 0.0461796049 0.593380420 0.36043997
## 178 0.4237399822 0.540480805 0.03577921
## 179 0.0002098446 0.007423691 0.99236646
## 180 0.2289460980 0.686893791 0.08416011
## 181 0.2289460980 0.686893791 0.08416011
## 182 0.4237399822 0.540480805 0.03577921
## 183 0.0078325591 0.216566804 0.77560064
## 184 0.0005195121 0.018173995 0.98130649
## 185 0.1070624427 0.707553600 0.18538396
## 186 0.0012855683 0.043764868 0.95494956
## 187 0.0031776352 0.101430219 0.89539215
## 188 0.4237399822 0.540480805 0.03577921
## 189 0.0078325591 0.216566804 0.77560064
## 190 0.0031776352 0.101430219 0.89539215
## 191 0.0031776352 0.101430219 0.89539215
## 192 0.0002098446 0.007423691 0.99236646
## 193 0.1070624427 0.707553600 0.18538396
## 194 0.0005195121 0.018173995 0.98130649
## 195 0.0031776352 0.101430219 0.89539215
## 196 0.0078325591 0.216566804 0.77560064
## 197 0.0191753167 0.398243691 0.58258099
## 198 0.0012855683 0.043764868 0.95494956
## 199 0.2289460980 0.686893791 0.08416011
## 200 0.0078325591 0.216566804 0.77560064

Asumsi Paralelisme dalam Regresi Logistik Ordinal

Asumsi paralelisme menyatakan bahwa efek dari variabel independen (misalnya, Waiting_Time dalam kasus Anda) terhadap log-odds adalah sama (konsisten) untuk semua ambang batas (cutpoints) antara kategori-kategori ordinal pada variabel dependen (misalnya, “Tidak Puas” ke “Cukup” dan “Cukup” ke “Puas”). Dengan kata lain, koefisien regresi untuk setiap variabel independen dianggap konstan di semua transisi antar-kategori.

Secara matematis, dalam model regresi logistik ordinal (proportional odds model), hubungan antara variabel independen \(X\) (Waiting_Time) dan probabilitas kumulatif dijelaskan sebagai:

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

  • \(Y\): Variabel dependen ordinal (Satisfaction: “Tidak Puas”, “Cukup”, “Puas”).
  • \(j\): Kategori ambang (misalnya, “Tidak Puas|Cukup” atau “Cukup|Puas”).
  • \(\alpha_j\): Intersep atau ambang batas untuk kategori \(j\).
  • \(\beta\): Koefisien regresi untuk variabel independen \(X\) (Waiting_Time), yang sama untuk semua ambang \(j\).

Pengujian Asumsi

# Membuat model logistik ordinal terlebih dahulu (menggunakan paket 'MASS')
library(MASS)
library(brant)
model <- polr(Satisfaction ~ Waiting_Time, data = df, Hess = TRUE)

# Uji Asumsi Proportional Odds
brant(model)
## -------------------------------------------- 
## Test for X2  df  probability 
## -------------------------------------------- 
## Omnibus      0.06    1   0.81
## Waiting_Time 0.06    1   0.81
## -------------------------------------------- 
## 
## H0: Parallel Regression Assumption holds

Kesimpulan : Diperoleh p-value (0.81) > 0.05 maka asumsi terpenuhi

Log Linear Model

Log-linear model adalah model statistik yang memodelkan frekuensi yang diharapkan (expected frequencies) dalam tabel kontingensi menggunakan fungsi logaritma. Model ini mengasumsikan bahwa logaritma dari frekuensi yang diharapkan dapat dinyatakan sebagai kombinasi linier dari efek utama (main effects) dan interaksi antar variabel kategorik.

Secara umum, log-linear model digunakan untuk:

Dalam log-linear model, frekuensi yang diharapkan (\(m_{ijk\ldots}\)) dalam sel tabel kontingensi dimodelkan sebagai:

\[ \log(m_{ijk\ldots}) = \mu + \lambda_i^A + \lambda_j^B + \lambda_k^C + \lambda_{ij}^{AB} + \lambda_{ik}^{AC} + \lambda_{jk}^{BC} + \lambda_{ijk}^{ABC} + \dots \]

2. Independensi dan Model Hierarkis

Log-linear model sering digunakan untuk menguji berbagai jenis independensi dalam tabel kontingensi:

Independensi Bersama (Mutual Independence)

Untuk tabel tiga arah \(X, Y, Z\), independensi bersama berarti: \[ P(X=i, Y=j, Z=k) = P(X=i)P(Y=j)P(Z=k) \] Model log-linear yang sesuai adalah: \[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z \] Tidak ada istilah interaksi dalam model ini, menunjukkan tidak ada asosiasi antar variabel.

Independensi Kondisional (Conditional Independence)

Misalnya, \(X\) dan \(Y\) independen diberikan \(Z\), ditulis: \[ P(X=i, Y=j | Z=k) = P(X=i | Z=k)P(Y=j | Z=k) \] Model log-linear yang sesuai adalah: \[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \] Tidak ada istilah \(\lambda_{ij}^{XY}\), menunjukkan tidak ada asosiasi langsung antara \(X\) dan \(Y\).

Model Hierarkis

Log-linear model biasanya hierarkis, artinya jika istilah interaksi tingkat tinggi disertakan (misalnya, \(\lambda_{ij}^{XY}\)), maka semua istilah tingkat rendah yang terkait (misalnya, \(\lambda_i^X, \lambda_j^Y\)) juga harus disertakan. Contoh model hierarkis untuk tabel tiga arah: - Model saturasi (memasukkan semua efek): \[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} + \lambda_{ijk}^{XYZ} \] - Model lebih sederhana (misalnya, tanpa interaksi tiga arah): \[ \log(\mu_{ijk}) = \lambda + \lambda_i^X + \lambda_j^Y + \lambda_k^Z + \lambda_{ij}^{XY} + \lambda_{ik}^{XZ} + \lambda_{jk}^{YZ} \]

Contoh
Dalam analisis data hukuman mati, model independensi bersama ditolak karena ada asosiasi signifikan. Model dengan independensi kondisional diuji untuk kecocokan.


3. Hubungan dengan Rasio Odds

Parameter interaksi dalam log-linear model memiliki hubungan langsung dengan rasio odds, yang juga dibahas dalam konteks tabel kontingensi dan regresi logistik:

Rasio Odds dalam Tabel 2x2

Untuk tabel 2x2, istilah interaksi \(\lambda_{ij}^{XY}\) terkait dengan log rasio odds: \[ \log(\theta) = \log\left(\frac{\mu_{11}\mu_{22}}{\mu_{12}\mu_{21}}\right) = \lambda_{11}^{XY} - \lambda_{12}^{XY} - \lambda_{21}^{XY} + \lambda_{22}^{XY} \] Ini memungkinkan interpretasi asosiasi dalam istilah rasio odds, yang lebih intuitif.

Tabel Multi-Dimensi

Dalam tabel tiga arah, rasio odds kondisional (misalnya, antara \(X\) dan \(Y\) diberikan \(Z\)) dapat dihitung dari istilah interaksi yang sesuai. Model tanpa istilah \(\lambda_{ij}^{XY}\) menyiratkan rasio odds kondisional sama dengan 1 (independensi kondisional).

Contoh
Dalam data hukuman mati, rasio odds kondisional antara hukuman dan ras terdakwa dihitung untuk setiap level ras korban, menunjukkan asosiasi yang signifikan.


4. Inferensi untuk Log-Linear Model

Inferensi statistik dalam log-linear model melibatkan estimasi parameter, pengujian kecocokan model, dan perbandingan model:

Estimasi Parameter

Parameter \(\lambda\) diestimasi menggunakan Maximum Likelihood Estimation (MLE) dengan asumsi distribusi Poisson. Perangkat lunak seperti SAS (PROC GENMOD) digunakan untuk perhitungan.

Goodness-of-Fit

Kecocokan model diuji menggunakan dua statistik utama: - Statistik Pearson Chi-Square: \[ X^2 = \sum \frac{(n_{ijk} - \hat{\mu}_{ijk})^2}{\hat{\mu}_{ijk}} \] - Statistik Rasio Likelihud (Deviansi): \[ G^2 = 2 \sum n_{ijk} \log\left(\frac{n_{ijk}}{\hat{\mu}_{ijk}}\right) \] Kedua statistik ini mengikuti distribusi chi-square dengan derajat kebebasan yang tergantung pada jumlah parameter dalam model. Nilai \(G^2\) atau \(X^2\) yang kecil menunjukkan kecocokan model yang baik.

Perbandingan Model

Model hierarkis dibandingkan menggunakan uji rasio likelihud untuk menentukan apakah istilah interaksi tertentu diperlukan. Misalnya, untuk menguji apakah interaksi \(\lambda_{ijk}^{XYZ} = 0\): \[ \Delta G^2 = G^2(\text{model sederhana}) - G^2(\text{model lengkap}) \] diuji terhadap distribusi chi-square dengan derajat kebebasan sama dengan selisih jumlah parameter.

Residu

Residu Pearson \(\frac{n_{ijk} - \hat{\mu}_{ijk}}{\sqrt{\hat{\mu}_{ijk}}}\) digunakan untuk mengidentifikasi sel yang tidak cocok dengan model, mirip dengan residu dalam uji chi-square.

Contoh
Dalam data kecelakaan mobil (tabel tiga arah: jenis kelamin, lokasi, penggunaan sabuk pengaman), model dengan independensi kondisional diuji dan dibandingkan dengan model saturasi. Deviansi menunjukkan model mana yang lebih cocok.


5. Model untuk Data Ordinal

Jika variabel dalam tabel kontingensi berskala ordinal, log-linear model dapat dimodifikasi untuk memanfaatkan urutan kategori:

Model Asosiasi Linear-by-Linear

Untuk tabel dua arah dengan variabel ordinal, model memasukkan istilah asosiasi linier: \[ \log(\mu_{ij}) = \lambda + \lambda_i^X + \lambda_j^Y + \beta u_i v_j \] di mana \(u_i\) dan \(v_j\) adalah skor ordinal untuk kategori \(X\) dan \(Y\), dan \(\beta\) mengukur kekuatan asosiasi linier. Model ini lebih sederhana daripada model saturasi dan memiliki daya uji lebih besar untuk tren linier.

Contoh
Dalam tabel yang menghubungkan tingkat pendidikan (ordinal) dan kepuasan kerja (ordinal), model linear-by-linear menunjukkan hubungan positif antara pendidikan yang lebih tinggi dan kepuasan yang lebih tinggi.

Ekstensi ke Tabel Multi-Dimensi

Untuk tabel tiga arah, model ordinal dapat mengasumsikan asosiasi homogen atau tren linier di setiap level variabel ketiga.


6. Model untuk Frekuensi Tingkat (Rate Data)

Log-linear model juga dapat digunakan untuk menganalisis data frekuensi tingkat (rate data), di mana frekuensi dihitung relatif terhadap jumlah paparan (exposure):

Model Poisson untuk Tingkat

\[ \log(\mu_{ij}) = \log(n_{ij}) + \lambda + \lambda_i^X + \lambda_j^Y + \lambda_{ij}^{XY} \] di mana \(\log(n_{ij})\) adalah offset untuk menyesuaikan jumlah paparan.

Contoh
Dalam analisis tingkat kematian akibat kanker paru-paru berdasarkan usia dan kebiasaan merokok, log-linear model digunakan untuk memodelkan tingkat kematian per 100.000 orang.


7. Estimasi dengan Sampel Kecil dan Data Jarang

Untuk tabel dengan frekuensi rendah atau sel kosong, estimasi MLE dalam log-linear model dapat bermasalah:

Data Jarang (Sparse Data)

Jika banyak sel memiliki frekuensi nol, estimasi parameter mungkin tidak stabil atau tak terbatas. Pendekatan seperti penalti atau metode Bayesian dapat digunakan.

Metode Iteratif

Algoritma seperti iterative proportional fitting (IPF) digunakan untuk mengestimasi frekuensi yang diharapkan dalam model hierarkis.

Contoh
Dalam tabel dengan banyak sel kosong (misalnya, data survei langka), model sederhana seperti independensi kondisional diuji untuk menghindari estimasi yang tidak stabil.

table_data <- matrix(c(40, 10, 30, 20), nrow=2,
  dimnames = list(Gender = c("Laki-Laki", "Perempuan"),
                  Gadget = c("Menyukai", "Tidak Menyukai")))
table_data
##            Gadget
## Gender      Menyukai Tidak Menyukai
##   Laki-Laki       40             30
##   Perempuan       10             20

Model Log Linear

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

Model Saturated

model_saturated <- loglm(~ Gender * Gadget, data = table_data)
summary(model_saturated)
## Formula:
## ~Gender * Gadget
## attr(,"variables")
## list(Gender, Gadget)
## attr(,"factors")
##        Gender Gadget Gender:Gadget
## Gender      1      0             1
## Gadget      0      1             1
## attr(,"term.labels")
## [1] "Gender"        "Gadget"        "Gender:Gadget"
## attr(,"order")
## [1] 1 1 2
## attr(,"intercept")
## [1] 1
## attr(,"response")
## [1] 0
## attr(,".Environment")
## <environment: R_GlobalEnv>
## 
## Statistics:
##                  X^2 df P(> X^2)
## Likelihood Ratio   0  0        1
## Pearson            0  0        1

Referensi

  1. https://bookdown.org/ssjackson300/ASM_Lecture_Notes/twocontingencytables.html
  2. https://dzchilds.github.io/stats-for-bio/contingency-tables.html
  3. Hosmer, D. W., & Lemeshow, S. (2013). Applied Logistic Regression (3rd ed.). Wiley.
  4. Kleinbaum, D. G., & Klein, M. (2010). Logistic Regression: A Self-Learning Text (3rd ed.). Springer.
  5. Roback, P., & Legler, J. (2021). Beyond Multiple Linear Regression: Applied Generalized Linear Models and Multilevel Models in R. [https://bookdown.org/roback/bookdown-BeyondMLR/ch-poissonreg.html]
  6. Wikipedia. (2023). Poisson Regression. [https://en.wikipedia.org/wiki/Poisson_regression]
  7. Frost, J. (2023). Poisson Regression Analysis: Overview & Example. [https://statisticsbyjim.com/regression/poisson-regression-analysis/]
  8. Generalized Linear Models Detailed Explanation https://dept.stat.lsa.umich.edu/~kshedden/Courses/Stat504/posts/glm/
  9. Generalized Linear Models Standard Textbook https://www.routledge.com/Generalized-Linear-Models/McCullagh-Nelder/p/book/9780412317606
  10. The Analysis Factor - Logistic Regression Models for Multinomial and Ordinal Variables
  11. UCLA Institute for Digital Research and Education - Ordinal Logistic Regression