1 Pendahuluan

Laporan ini menganalisis berbagai aspek prediktif perilaku pelanggan, kinerja akademik, perputaran karyawan, dan kunjungan layanan kesehatan menggunakan empat dataset sekunder yang berbeda. Setiap dataset diaplikasikan pada satu jenis model statistik yang spesifik karena struktur dan bentuk variabel responnya yang berbeda, mulai dari data biner, ordinal, multinomial, hingga data diskret (cacah).

Empat model yang digunakan adalah:

  1. Regresi logistik biner Digunakan untuk menganalisis dan memprediksi loyalitas pelanggan telekomunikasi. Variabel responnya bersifat biner (dikotomi), yaitu apakah pelanggan berhenti berlangganan (Churn = Yes) atau bertahan (Churn = No), berdasarkan profil layanan dan biaya bulanan mereka.
  2. Regresi logistik ordinal Digunakan untuk menganalisis performa akademik siswa yang diukur melalui nilai akhir (G3). Variabel respon ini dikelompokkan ke dalam skala ordinal (Rendah < Sedang < Tinggi) untuk melihat pengaruh latar belakang keluarga, konsumsi alkohol, dan waktu belajar terhadap pencapaian siswa.
  3. Regresi logistik multinomial Digunakan untuk menganalisis peran jabatan karyawan (JobRole) di perusahaan. Variabel respon memiliki banyak kategori tanpa urutan, sehingga setiap kategori peran akan dibandingkan terhadap satu kategori acuan (baseline).
  4. Regresi Poisson Digunakan untuk menganalisis data cacah berupa jumlah dokter yang dikunjungi oleh pasien lansia (Number of Doctors Visited). Model ini sangat cocok karena datanya berupa bilangan bulat non-negatif (\(0, 1, 2, \dots\)) yang dipengaruhi oleh faktor kesehatan fisik, mental, dan gangguan tidur pasien.

1.1 Istilah Penting dalam Bahasa Sederhana

Istilah Penjelasan dalam laporan ini
Churn / Retensi Kondisi ketika pelanggan memutuskan untuk berhenti berlangganan layanan (Telco Customer Churn). Dalam data biner, 1 (Yes) berarti churn (berhenti) dan 0 (No) berarti loyal (bertahan).
Performa G3 Nilai akhir siswa pada mata pelajaran bahasa Portugis (student-por.csv). Dalam model ordinal, nilai ini dikelompokkan menjadi tiga tingkatan: Rendah, Sedang, dan Tinggi.
JobRole Peran jabatan karyawan dalam dataset HR Employee Attrition. Digunakan sebagai respon multinomial karena memiliki banyak kategori tanpa urutan.
Data Cacah (Count Data) Data yang nilainya berupa bilangan bulat non-negatif hasil dari penghitungan. Dalam laporan ini, contohnya adalah jumlah dokter yang dikunjungi oleh pasien lansia.
Odds Ratio (OR) Ukuran untuk membaca pengaruh variabel prediktor pada model logistik (biner dan ordinal). Untuk model ordinal dengan fungsi MASS::polr di R, nilai \(\text{OR} = \exp(\beta)\) diinterpretasikan sebagai perubahan peluang (odds) untuk berada pada kategori performa/risiko yang lebih tinggi, karena spesifikasi modelnya menggunakan formula:\[\text{logit}(P(Y \le j)) = \alpha_j - x\beta\]
Relative risk ratio / RRR Ukuran untuk membaca pengaruh prediktor pada regresi multinomial. RRR membandingkan satu kategori respon terhadap kategori acuan.
Incidence rate ratio / IRR Ukuran untuk membaca pengaruh prediktor pada regresi Poisson. IRR lebih dari 1 berarti rata-rata jumlah kejadian lebih tinggi.
Data leakage Kondisi ketika model diberi variabel yang sebenarnya sudah mengandung jawaban. Variabel seperti ini harus dihindari agar model tidak terlihat bagus secara palsu.

2 Sumber Data dan Rancangan Analisis

Dataset yang digunakan bersumber dari repositori publik/sekunder: IBM Business Analytics Community untuk data perputaran pelanggan telekomunikasi (Telco Customer Churn) dan perputaran karyawan (HR Employee Attrition), UCI Machine Learning Repository untuk data performa akademik siswa bahasa Portugis (student-por), serta National Poll on Healthy Aging (NPHA) untuk data kunjungan dokter pasien lansia.

Model Dataset Variabel respon Bentuk respon Tujuan analisis
Logistik biner Telco Customer Churn Churn 0/1 Menjelaskan peluang pelanggan berhenti berlangganan
Ordinal Student Performance (Portugis) Performa_G3 Rendah < Sedang < Tinggi Menjelaskan faktor penentu tingkat capaian akademik siswa
Multinomial HR Employee Attrition JobRole Banyak kategori tanpa urutan Menjelaskan kecenderungan penempatan atau peran jabatan karyawan
Poisson NPHA Doctor Visits Number of Doctors Visited Hitungan/count Menjelaskan frekuensi jumlah kunjungan pasien ke dokter

3 Konsep Singkat Model

3.1 Regresi Logistik Biner

Regresi logistik biner dipakai ketika variabel respon hanya memiliki dua kategori (dikotomi). Dalam laporan ini, responnya adalah apakah pelanggan berhenti berlangganan atau tidak (Churn = Yes atau Churn = No).

\[ \log\left(\frac{P(Y=1)}{P(Y=0)}\right) = \beta_0 + \beta_1X_1 + \cdots + \beta_pX_p \]

3.2 Regresi Logistik Ordinal

Regresi ordinal dipakai ketika variabel respon memiliki tingkatan atau urutan yang jelas. Dalam laporan ini, responnya adalah performa nilai akhir siswa (G3) yang dikelompokkan menjadi kategori Rendah < Sedang < Tinggi.

Pada analisis ini model diestimasi menggunakan fungsi MASS::polr. Spesifikasi model polr adalah:

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

Artinya, tanda koefisien harus dibaca dengan hati-hati. Jika koefisien \(\beta\) bernilai positif, maka nilai \(\alpha_j - x\beta\) menjadi lebih kecil. Akibatnya, peluang kumulatif untuk berada pada kategori rendah \(P(Y \le j)\) menurun, sehingga peluang untuk berada pada kategori yang lebih tinggi meningkat. Dengan urutan tingkat nilai Rendah < Sedang < Tinggi, koefisien positif berarti variabel prediktor tersebut mendorong kecenderungan siswa untuk bergeser ke tingkat performa akademik yang lebih tinggi/baik.

3.3 Regresi Logistik Multinomial

Regresi multinomial dipakai ketika variabel respon memiliki lebih dari dua kategori dan kategorinya tidak memiliki urutan atau tingkatan (nominal). Dalam laporan ini, respon yang digunakan adalah jenis peran jabatan karyawan (JobRole), di mana setiap kategori peran akan dibandingkan terhadap satu kategori acuan (baseline / reference cell).

3.4 Regresi Poisson

Regresi Poisson dipakai ketika variabel respon berupa data cacah atau jumlah frekuensi kejadian (count data). Dalam laporan ini, responnya adalah jumlah dokter yang dikunjungi oleh pasien lansia. Karena data berbentuk individual-level (satu baris per pasien), model Poisson digunakan tanpa offset.

\[ Y_i \sim Poisson(\mu_i) \]

\[ \log(\mu_i) = \beta_0 + \beta_1X_{1i} + \cdots + \beta_pX_{pi} \]

Berbeda dengan data kelompok/agregat, data tingkat individu ini menganalisis rata-rata nilai harapan (\(\mu_i\)) jumlah kunjungan dokter secara langsung berdasarkan kondisi kesehatan fisik, mental, dan gangguan tidur pasien tanpa memerlukan variabel offset.

4 Persiapan Package

required_packages <- c(
  "readxl", "readr", "dplyr", "tidyr", "ggplot2",
  "MASS", "nnet", "broom", "knitr", "scales", "stringr", "forcats", "tibble", "purrr",
  "pROC", "ordinal"
)

missing_packages <- required_packages[!(required_packages %in% rownames(installed.packages()))]
if (length(missing_packages) > 0) {
  install.packages(missing_packages, repos = "https://cloud.r-project.org")
}

invisible(lapply(required_packages, library, character.only = TRUE))

fmt_num <- function(x, digits = 3) {
  ifelse(is.na(x), NA_character_, format(round(x, digits), nsmall = digits, big.mark = ","))
}

fmt_p <- function(x) {
  ifelse(is.na(x), NA_character_, ifelse(x < 0.001, "< 0.001", fmt_num(x, 3)))
}

ratio_phrase <- function(ratio, measure, up_text, down_text) {
  if (is.na(ratio)) return("nilai rasio tidak tersedia")
  if (ratio > 1) {
    paste0(measure, " = ", fmt_num(ratio, 3), ", sehingga ", up_text)
  } else if (ratio < 1) {
    paste0(measure, " = ", fmt_num(ratio, 3), ", sehingga ", down_text)
  } else {
    paste0(measure, " = 1, sehingga tidak menunjukkan perubahan kecenderungan")
  }
}

modus <- function(x) {
  ux <- na.omit(unique(x))
  ux[which.max(tabulate(match(x, ux)))]
}

# Label untuk model biner (Telco Customer Churn)
label_biner <- function(term) {
  dplyr::case_when(
    term == "tenure"              ~ "durasi berlangganan (bulan)",
    term == "MonthlyCharges"      ~ "biaya bulanan pelanggan",
    term == "TotalCharges"        ~ "total biaya kumulatif pelanggan",
    term == "SeniorCitizen"       ~ "status lansia (1 = Ya, 0 = Tidak)",

    stringr::str_starts(term, "gender")         ~ paste0("jenis kelamin: ", stringr::str_remove(term, "gender"), " dibanding kategori acuan"),
    stringr::str_starts(term, "Contract")       ~ paste0("jenis kontrak: ", stringr::str_remove(term, "Contract"), " dibanding kategori acuan"),
    stringr::str_starts(term, "InternetService")~ paste0("layanan internet: ", stringr::str_remove(term, "InternetService"), " dibanding kategori acuan"),
    stringr::str_starts(term, "PaymentMethod")  ~ paste0("metode pembayaran: ", stringr::str_remove(term, "PaymentMethod"), " dibanding kategori acuan"),
    stringr::str_starts(term, "OnlineSecurity") ~ paste0("keamanan online: ", stringr::str_remove(term, "OnlineSecurity"), " dibanding kategori acuan"),
    stringr::str_starts(term, "OnlineBackup")   ~ paste0("pencadangan online: ", stringr::str_remove(term, "OnlineBackup"), " dibanding kategori acuan"),
    stringr::str_starts(term, "DeviceProtection")~ paste0("perlindungan perangkat: ", stringr::str_remove(term, "DeviceProtection"), " dibanding kategori acuan"),
    stringr::str_starts(term, "TechSupport")    ~ paste0("dukungan teknis: ", stringr::str_remove(term, "TechSupport"), " dibanding kategori acuan"),
    stringr::str_starts(term, "StreamingTV")    ~ paste0("layanan streaming TV: ", stringr::str_remove(term, "StreamingTV"), " dibanding kategori acuan"),
    stringr::str_starts(term, "StreamingMovies")~ paste0("layanan streaming film: ", stringr::str_remove(term, "StreamingMovies"), " dibanding kategori acuan"),
    stringr::str_starts(term, "MultipleLines")  ~ paste0("multiple lines: ", stringr::str_remove(term, "MultipleLines"), " dibanding kategori acuan"),
    stringr::str_starts(term, "PhoneService")   ~ paste0("layanan telepon: ", stringr::str_remove(term, "PhoneService"), " dibanding kategori acuan"),
    stringr::str_starts(term, "Partner")        ~ paste0("status pasangan: ", stringr::str_remove(term, "Partner"), " dibanding kategori acuan"),
    stringr::str_starts(term, "Dependents")     ~ paste0("status tanggungan: ", stringr::str_remove(term, "Dependents"), " dibanding kategori acuan"),
    stringr::str_starts(term, "PaperlessBilling")~ paste0("tagihan tanpa kertas: ", stringr::str_remove(term, "PaperlessBilling"), " dibanding kategori acuan"),
    TRUE ~ term
  )
}

# Label untuk model ordinal (Student Performance)
label_ordinal <- function(term) {
  dplyr::case_when(
    term == "studytime"  ~ "waktu belajar mingguan (skala 1-4)",
    term == "failures"   ~ "jumlah kegagalan kelas sebelumnya",
    term == "Dalc"       ~ "konsumsi alkohol hari kerja (skala 1-5)",
    term == "Walc"       ~ "konsumsi alkohol akhir pekan (skala 1-5)",
    term == "age"        ~ "usia siswa",
    term == "absences"   ~ "jumlah absensi/ketidakhadiran sekolah",
    term == "traveltime" ~ "waktu perjalanan ke sekolah",
    term == "famrel"     ~ "kualitas hubungan keluarga",
    term == "freetime"   ~ "waktu luang setelah sekolah",
    term == "goout"      ~ "frekuensi pergi bermain bersama teman",
    term == "health"     ~ "status kesehatan saat ini",
    term == "Medu"       ~ "tingkat pendidikan ibu",
    term == "Fedu"       ~ "tingkat pendidikan ayah",

    stringr::str_starts(term, "school")   ~ paste0("sekolah: ", stringr::str_remove(term, "school"), " dibanding kategori acuan"),
    stringr::str_starts(term, "sex")      ~ paste0("jenis kelamin siswa: ", stringr::str_remove(term, "sex"), " dibanding kategori acuan"),
    stringr::str_starts(term, "address")  ~ paste0("tipe alamat tinggal: ", stringr::str_remove(term, "address"), " dibanding kategori acuan"),
    stringr::str_starts(term, "Mjob")     ~ paste0("pekerjaan ibu: ", stringr::str_remove(term, "Mjob"), " dibanding kategori acuan"),
    stringr::str_starts(term, "Fjob")     ~ paste0("pekerjaan ayah: ", stringr::str_remove(term, "Fjob"), " dibanding kategori acuan"),
    stringr::str_starts(term, "internet") ~ paste0("akses internet di rumah: ", stringr::str_remove(term, "internet"), " dibanding acuan"),
    stringr::str_starts(term, "romantic") ~ paste0("memiliki hubungan romantis: ", stringr::str_remove(term, "romantic"), " dibanding acuan"),
    stringr::str_starts(term, "famsup")   ~ paste0("dukungan pendidikan keluarga: ", stringr::str_remove(term, "famsup"), " dibanding acuan"),
    stringr::str_starts(term, "schoolsup")~ paste0("dukungan pendidikan sekolah: ", stringr::str_remove(term, "schoolsup"), " dibanding acuan"),
    stringr::str_starts(term, "paid")     ~ paste0("les berbayar: ", stringr::str_remove(term, "paid"), " dibanding acuan"),
    stringr::str_starts(term, "activities")~paste0("kegiatan ekstrakulikuler: ", stringr::str_remove(term, "activities"), " dibanding acuan"),
    stringr::str_starts(term, "higher")   ~ paste0("ingin melanjutkan pendidikan tinggi: ", stringr::str_remove(term, "higher"), " dibanding acuan"),
    stringr::str_starts(term, "nursery")  ~ paste0("pernah ikut taman kanak-kanak: ", stringr::str_remove(term, "nursery"), " dibanding acuan"),
    stringr::str_starts(term, "famsize")  ~ paste0("ukuran keluarga: ", stringr::str_remove(term, "famsize"), " dibanding acuan"),
    stringr::str_starts(term, "Pstatus")  ~ paste0("status orang tua: ", stringr::str_remove(term, "Pstatus"), " dibanding acuan"),
    stringr::str_starts(term, "reason")   ~ paste0("alasan pilih sekolah: ", stringr::str_remove(term, "reason"), " dibanding acuan"),
    stringr::str_starts(term, "guardian") ~ paste0("wali siswa: ", stringr::str_remove(term, "guardian"), " dibanding acuan"),
    TRUE ~ term
  )
}

# Label untuk model multinomial (HR Employee Attrition - JobRole)
label_multinom <- function(term) {
  dplyr::case_when(
    term == "Age"                    ~ "usia karyawan",
    term == "MonthlyIncome"          ~ "pendapatan bulanan karyawan",
    term == "DistanceFromHome"       ~ "jarak dari rumah ke kantor (km)",
    term == "PercentSalaryHike"      ~ "persentase kenaikan gaji",
    term == "TotalWorkingYears"      ~ "total masa kerja (tahun)",
    term == "YearsAtCompany"         ~ "lama bekerja di perusahaan saat ini",
    term == "YearsInCurrentRole"     ~ "lama memegang peran jabatan saat ini",
    term == "JobSatisfaction"        ~ "tingkat kepuasan kerja (skala 1-4)",
    term == "WorkLifeBalance"        ~ "keseimbangan hidup-kerja (skala 1-4)",
    term == "EnvironmentSatisfaction"~ "kepuasan lingkungan kerja (skala 1-4)",
    term == "JobInvolvement"         ~ "keterlibatan dalam pekerjaan (skala 1-4)",
    term == "JobLevel"               ~ "level jabatan karyawan",
    term == "NumCompaniesWorked"     ~ "jumlah perusahaan sebelumnya",
    term == "TrainingTimesLastYear"  ~ "frekuensi pelatihan tahun lalu",
    term == "StockOptionLevel"       ~ "level opsi saham karyawan",
    term == "YearsSinceLastPromotion"~ "tahun sejak promosi terakhir",
    term == "YearsWithCurrManager"   ~ "lama dengan manajer saat ini",
    term == "DailyRate"              ~ "tarif harian karyawan",
    term == "HourlyRate"             ~ "tarif per jam karyawan",
    term == "MonthlyRate"            ~ "tarif bulanan karyawan",

    stringr::str_starts(term, "Attrition")           ~ paste0("status keluar karyawan: ", stringr::str_remove(term, "Attrition"), " dibanding acuan"),
    stringr::str_starts(term, "BusinessTravel")      ~ paste0("frekuensi perjalanan bisnis: ", stringr::str_remove(term, "BusinessTravel"), " dibanding acuan"),
    stringr::str_starts(term, "Department")          ~ paste0("departemen kerja: ", stringr::str_remove(term, "Department"), " dibanding acuan"),
    stringr::str_starts(term, "EducationField")      ~ paste0("bidang pendidikan: ", stringr::str_remove(term, "EducationField"), " dibanding acuan"),
    stringr::str_starts(term, "Gender")              ~ paste0("jenis kelamin karyawan: ", stringr::str_remove(term, "Gender"), " dibanding acuan"),
    stringr::str_starts(term, "MaritalStatus")       ~ paste0("status pernikahan: ", stringr::str_remove(term, "MaritalStatus"), " dibanding acuan"),
    stringr::str_starts(term, "OverTime")            ~ paste0("bekerja lembur: ", stringr::str_remove(term, "OverTime"), " dibanding acuan"),
    TRUE ~ term
  )
}

# Label untuk model Poisson (NPHA Doctor Visits)
label_poisson <- function(term) {
  dplyr::case_when(
    term == "Age"                                    ~ "kelompok usia pasien",
    term == "Phyiscal.Health"                        ~ "status kesehatan fisik pasien",
    term == "Mental.Health"                          ~ "status kesehatan mental pasien",
    term == "Dental.Health"                          ~ "status kesehatan gigi pasien",
    term == "Employment"                             ~ "status pekerjaan pasien",
    term == "Trouble.Sleeping"                       ~ "tingkat kesulitan tidur umum",
    term == "Stress.Keeps.Patient.from.Sleeping"     ~ "gangguan tidur akibat stres",
    term == "Medication.Keeps.Patient.from.Sleeping" ~ "gangguan tidur akibat efek obat",
    term == "Pain.Keeps.Patient.from.Sleeping"       ~ "gangguan tidur akibat rasa nyeri fisik",
    term == "Bathroom.Needs.Keeps.Patient.from.Sleeping" ~ "gangguan tidur akibat kebutuhan ke toilet",
    term == "Uknown.Keeps.Patient.from.Sleeping"     ~ "gangguan tidur akibat penyebab tidak diketahui",
    term == "Prescription.Sleep.Medication"          ~ "penggunaan obat tidur resep dokter",
    term == "Race"                                   ~ "kategori ras pasien",
    term == "Gender"                                 ~ "jenis kelamin pasien",
    TRUE ~ term
  )
}

5 Import Data

data_dir <- "."
file_biner   <- file.path(data_dir, "WA_Fn-UseC_-Telco-Customer-Churn.csv")
file_ordinal <- file.path(data_dir, "student-por.csv")
file_multinom<- file.path(data_dir, "WA_Fn-UseC_-HR-Employee-Attrition.csv")
file_poisson <- file.path(data_dir, "NPHA-doctor-visits.csv")

check_file <- function(path) {
  if (!file.exists(path)) {
    stop(paste("File tidak ditemukan:", path,
               "\nPastikan file Rmd dan dataset berada di folder yang sama atau ubah data_dir."))
  }
}

invisible(lapply(c(file_biner, file_ordinal, file_multinom, file_poisson), check_file))

data_biner_raw   <- readr::read_csv(file_biner,   show_col_types = FALSE)
data_ordinal_raw <- readr::read_csv(file_ordinal,  show_col_types = FALSE)
data_multinom_raw<- readr::read_csv(file_multinom, show_col_types = FALSE)
data_poisson_raw <- readr::read_csv(file_poisson,  show_col_types = FALSE)
dimensi_data <- tibble::tibble(
  Dataset = c("Telco Customer Churn", "Student Performance (Portugis)",
              "HR Employee Attrition", "NPHA Doctor Visits"),
  Jumlah_Baris = c(nrow(data_biner_raw), nrow(data_ordinal_raw),
                   nrow(data_multinom_raw), nrow(data_poisson_raw)),
  Jumlah_Kolom = c(ncol(data_biner_raw), ncol(data_ordinal_raw),
                   ncol(data_multinom_raw), ncol(data_poisson_raw))
)
knitr::kable(dimensi_data, caption = "Dimensi Data Awal")
Dimensi Data Awal
Dataset Jumlah_Baris Jumlah_Kolom
Telco Customer Churn 7043 21
Student Performance (Portugis) 649 33
HR Employee Attrition 1470 35
NPHA Doctor Visits 714 15

6 Bagian I — Regresi Logistik Biner

6.1 Tujuan Analisis

Tujuan analisis ini adalah mengetahui faktor-faktor yang berhubungan dengan peluang pelanggan telekomunikasi melakukan churn (berhenti berlangganan). Churn berarti pelanggan memutuskan untuk tidak melanjutkan layanan dari perusahaan.

6.2 Variabel yang Digunakan

Variabel respon adalah Churn:

Variabel yang Dipakai pada Model Biner
Nilai Churn Arti sederhana
No (0) Pelanggan tidak berhenti berlangganan (loyal)
Yes (1) Pelanggan berhenti berlangganan (churn)
Variabel Peran Alasan_dipakai
Churn Respon Mewakili status churn pelanggan (biner: Yes/No)
gender Prediktor Karakteristik demografi pelanggan
SeniorCitizen Prediktor Status lansia dapat mempengaruhi loyalitas pelanggan
Partner Prediktor Pelanggan dengan pasangan mungkin lebih stabil
Dependents Prediktor Tanggungan keluarga dapat mempengaruhi keputusan churn
tenure Prediktor Semakin lama berlangganan, biasanya semakin loyal
PhoneService Prediktor Penggunaan layanan telepon
MultipleLines Prediktor Penggunaan lebih dari satu jalur telepon
InternetService Prediktor Jenis layanan internet yang digunakan
OnlineSecurity Prediktor Penggunaan layanan keamanan online
OnlineBackup Prediktor Penggunaan layanan pencadangan online
DeviceProtection Prediktor Penggunaan layanan perlindungan perangkat
TechSupport Prediktor Penggunaan layanan dukungan teknis
StreamingTV Prediktor Penggunaan layanan streaming TV
StreamingMovies Prediktor Penggunaan layanan streaming film
Contract Prediktor Jenis kontrak berpengaruh kuat terhadap churn
PaperlessBilling Prediktor Metode tagihan
PaymentMethod Prediktor Metode pembayaran yang digunakan
MonthlyCharges Prediktor Biaya bulanan yang dibayar pelanggan
TotalCharges Prediktor Total biaya kumulatif selama berlangganan

Variabel customerID tidak digunakan karena hanya merupakan identitas unik pelanggan dan tidak memiliki nilai prediktif. TotalCharges perlu diperiksa tipe datanya karena terkadang terbaca sebagai karakter akibat nilai kosong.

6.3 Cleaning Data Biner

data_biner <- data_biner_raw %>%
  dplyr::mutate(
    TotalCharges = as.numeric(TotalCharges),
    Churn_num    = ifelse(Churn == "Yes", 1L, 0L),
    Churn_label  = factor(Churn, levels = c("No", "Yes")),
    gender            = factor(gender),
    Partner           = factor(Partner),
    Dependents        = factor(Dependents),
    PhoneService      = factor(PhoneService),
    MultipleLines     = factor(MultipleLines),
    InternetService   = factor(InternetService),
    OnlineSecurity    = factor(OnlineSecurity),
    OnlineBackup      = factor(OnlineBackup),
    DeviceProtection  = factor(DeviceProtection),
    TechSupport       = factor(TechSupport),
    StreamingTV       = factor(StreamingTV),
    StreamingMovies   = factor(StreamingMovies),
    Contract          = factor(Contract),
    PaperlessBilling  = factor(PaperlessBilling),
    PaymentMethod     = factor(PaymentMethod)
  ) %>%
  tidyr::drop_na()

knitr::kable(head(data_biner %>%
  dplyr::select(customerID, gender, SeniorCitizen, tenure, Contract,
                MonthlyCharges, TotalCharges, Churn_label), 6),
  caption = "Cuplikan Data Biner Setelah Cleaning")
Cuplikan Data Biner Setelah Cleaning
customerID gender SeniorCitizen tenure Contract MonthlyCharges TotalCharges Churn_label
7590-VHVEG Female 0 1 Month-to-month 29.85 29.85 No
5575-GNVDE Male 0 34 One year 56.95 1889.50 No
3668-QPYBK Male 0 2 Month-to-month 53.85 108.15 Yes
7795-CFOCW Male 0 45 One year 42.30 1840.75 No
9237-HQITU Female 0 2 Month-to-month 70.70 151.65 Yes
9305-CDSKC Female 0 8 Month-to-month 99.65 820.50 Yes

Variabel TotalCharges dikonversi ke numerik. Baris dengan nilai kosong dihapus menggunakan drop_na(). Variabel kategorik seperti gender, Contract, InternetService, dan PaymentMethod diubah menjadi factor agar R dapat membuat dummy variable secara otomatis.

6.4 Analisis Deskriptif Data Biner

dist_biner <- data_biner %>%
  dplyr::count(Churn_label) %>%
  dplyr::mutate(Persen = scales::percent(n / sum(n)))
knitr::kable(dist_biner, caption = "Distribusi Status Churn Pelanggan")
Distribusi Status Churn Pelanggan
Churn_label n Persen
No 5163 73%
Yes 1869 27%
summary_biner <- data_biner %>%
  dplyr::summarise(
    Rata_rata_tenure       = mean(tenure),
    Median_tenure          = median(tenure),
    Rata_rata_MonthlyCharges = mean(MonthlyCharges),
    Median_MonthlyCharges  = median(MonthlyCharges)
  )
knitr::kable(summary_biner, caption = "Ringkasan Prediktor Utama Data Biner")
Ringkasan Prediktor Utama Data Biner
Rata_rata_tenure Median_tenure Rata_rata_MonthlyCharges Median_MonthlyCharges
32.42179 29 64.79821 70.35
ggplot(data_biner, aes(x = Churn_label)) +
  geom_bar(fill = c("#2f7f73", "#b84f5a")) +
  labs(
    title = "Distribusi Status Churn Pelanggan Telekomunikasi",
    x     = "Status pelanggan",
    y     = "Jumlah pelanggan"
  ) +
  theme_minimal()

6.5 Pembentukan Model Regresi Logistik Biner

Model biner dibentuk untuk menjelaskan peluang Churn = Yes, yaitu pelanggan berhenti berlangganan.

\[ \log\left(\frac{P(\text{Churn}=1)}{P(\text{Churn}=0)}\right) = \beta_0 + \beta_1\,\text{gender} + \beta_2\,\text{SeniorCitizen} + \cdots + \beta_p X_p \]

model_biner <- glm(
  Churn_num ~ gender + SeniorCitizen + Partner + Dependents + tenure +
    PhoneService + MultipleLines + InternetService +
    OnlineSecurity + OnlineBackup + DeviceProtection + TechSupport +
    StreamingTV + StreamingMovies + Contract + PaperlessBilling +
    PaymentMethod + MonthlyCharges + TotalCharges,
  data   = data_biner,
  family = binomial(link = "logit")
)

summary(model_biner)
## 
## Call:
## glm(formula = Churn_num ~ gender + SeniorCitizen + Partner + 
##     Dependents + tenure + PhoneService + MultipleLines + InternetService + 
##     OnlineSecurity + OnlineBackup + DeviceProtection + TechSupport + 
##     StreamingTV + StreamingMovies + Contract + PaperlessBilling + 
##     PaymentMethod + MonthlyCharges + TotalCharges, family = binomial(link = "logit"), 
##     data = data_biner)
## 
## Coefficients: (7 not defined because of singularities)
##                                         Estimate  Std. Error z value
## (Intercept)                           1.16528747  0.81513548   1.430
## genderMale                           -0.02183273  0.06480439  -0.337
## SeniorCitizen                         0.21677504  0.08453049   2.564
## PartnerYes                           -0.00038400  0.07782976  -0.005
## DependentsYes                        -0.14848782  0.08973115  -1.655
## tenure                               -0.06058758  0.00623568  -9.716
## PhoneServiceYes                       0.17146779  0.64869215   0.264
## MultipleLinesNo phone service                 NA          NA      NA
## MultipleLinesYes                      0.44839544  0.17725585   2.530
## InternetServiceFiber optic            1.74747491  0.79807951   2.190
## InternetServiceNo                    -1.78629472  0.80726809  -2.213
## OnlineSecurityNo internet service             NA          NA      NA
## OnlineSecurityYes                    -0.20542004  0.17868789  -1.150
## OnlineBackupNo internet service               NA          NA      NA
## OnlineBackupYes                       0.02604185  0.17540119   0.148
## DeviceProtectionNo internet service           NA          NA      NA
## DeviceProtectionYes                   0.14737500  0.17637397   0.836
## TechSupportNo internet service                NA          NA      NA
## TechSupportYes                       -0.18049681  0.18060228  -0.999
## StreamingTVNo internet service                NA          NA      NA
## StreamingTVYes                        0.59050741  0.32630948   1.810
## StreamingMoviesNo internet service            NA          NA      NA
## StreamingMoviesYes                    0.59929571  0.32668370   1.834
## ContractOne year                     -0.66079529  0.10758532  -6.142
## ContractTwo year                     -1.35710624  0.17644517  -7.691
## PaperlessBillingYes                   0.34235364  0.07449538   4.596
## PaymentMethodCredit card (automatic) -0.08779182  0.11407927  -0.770
## PaymentMethodElectronic check         0.30446727  0.09449652   3.222
## PaymentMethodMailed check            -0.05758719  0.11491136  -0.501
## MonthlyCharges                       -0.04034353  0.03175504  -1.270
## TotalCharges                          0.00032894  0.00007063   4.657
##                                                  Pr(>|z|)    
## (Intercept)                                       0.15284    
## genderMale                                        0.73619    
## SeniorCitizen                                     0.01033 *  
## PartnerYes                                        0.99606    
## DependentsYes                                     0.09796 .  
## tenure                               < 0.0000000000000002 ***
## PhoneServiceYes                                   0.79153    
## MultipleLinesNo phone service                          NA    
## MultipleLinesYes                                  0.01142 *  
## InternetServiceFiber optic                        0.02855 *  
## InternetServiceNo                                 0.02691 *  
## OnlineSecurityNo internet service                      NA    
## OnlineSecurityYes                                 0.25031    
## OnlineBackupNo internet service                        NA    
## OnlineBackupYes                                   0.88197    
## DeviceProtectionNo internet service                    NA    
## DeviceProtectionYes                               0.40339    
## TechSupportNo internet service                         NA    
## TechSupportYes                                    0.31759    
## StreamingTVNo internet service                         NA    
## StreamingTVYes                                    0.07035 .  
## StreamingMoviesNo internet service                     NA    
## StreamingMoviesYes                                0.06658 .  
## ContractOne year                       0.0000000008145911 ***
## ContractTwo year                       0.0000000000000146 ***
## PaperlessBillingYes                    0.0000043143131817 ***
## PaymentMethodCredit card (automatic)              0.44156    
## PaymentMethodElectronic check                     0.00127 ** 
## PaymentMethodMailed check                         0.61627    
## MonthlyCharges                                    0.20392    
## TotalCharges                           0.0000032026251598 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 8143.4  on 7031  degrees of freedom
## Residual deviance: 5826.3  on 7008  degrees of freedom
## AIC: 5874.3
## 
## Number of Fisher Scoring iterations: 6

6.6 Persamaan Model Hasil Estimasi Biner

Chunk model-biner di atas sudah membangun model regresi logistik biner menggunakan fungsi glm() dengan family = binomial(link = "logit"). Persamaan umum hasil estimasi adalah:

\[ \widehat{\eta_i} = \log\left(\frac{\widehat{P}(\text{Churn}_i=1)}{\widehat{P}(\text{Churn}_i=0)}\right) = \hat{\beta}_0 + \hat{\beta}_1 X_{1i} + \cdots + \hat{\beta}_p X_{pi} \]

Dengan:

\[ \widehat{P}(\text{Churn}_i=1) = \frac{\exp(\widehat{\eta_i})}{1+\exp(\widehat{\eta_i})} \]

Parameter yang Membentuk Persamaan Regresi Logistik Biner
Variabel_Mudah term estimate Odds_Ratio
(Intercept) (Intercept) 1.1653 3.2068
jenis kelamin: Male dibanding kategori acuan genderMale -0.0218 0.9784
status lansia (1 = Ya, 0 = Tidak) SeniorCitizen 0.2168 1.2421
status pasangan: Yes dibanding kategori acuan PartnerYes -0.0004 0.9996
status tanggungan: Yes dibanding kategori acuan DependentsYes -0.1485 0.8620
durasi berlangganan (bulan) tenure -0.0606 0.9412
layanan telepon: Yes dibanding kategori acuan PhoneServiceYes 0.1715 1.1870
multiple lines: No phone service dibanding kategori acuan MultipleLinesNo phone service NA NA
multiple lines: Yes dibanding kategori acuan MultipleLinesYes 0.4484 1.5658
layanan internet: Fiber optic dibanding kategori acuan InternetServiceFiber optic 1.7475 5.7401
layanan internet: No dibanding kategori acuan InternetServiceNo -1.7863 0.1676
keamanan online: No internet service dibanding kategori acuan OnlineSecurityNo internet service NA NA
keamanan online: Yes dibanding kategori acuan OnlineSecurityYes -0.2054 0.8143
pencadangan online: No internet service dibanding kategori acuan OnlineBackupNo internet service NA NA
pencadangan online: Yes dibanding kategori acuan OnlineBackupYes 0.0260 1.0264
perlindungan perangkat: No internet service dibanding kategori acuan DeviceProtectionNo internet service NA NA
perlindungan perangkat: Yes dibanding kategori acuan DeviceProtectionYes 0.1474 1.1588
dukungan teknis: No internet service dibanding kategori acuan TechSupportNo internet service NA NA
dukungan teknis: Yes dibanding kategori acuan TechSupportYes -0.1805 0.8349
layanan streaming TV: No internet service dibanding kategori acuan StreamingTVNo internet service NA NA
layanan streaming TV: Yes dibanding kategori acuan StreamingTVYes 0.5905 1.8049
layanan streaming film: No internet service dibanding kategori acuan StreamingMoviesNo internet service NA NA
layanan streaming film: Yes dibanding kategori acuan StreamingMoviesYes 0.5993 1.8208
jenis kontrak: One year dibanding kategori acuan ContractOne year -0.6608 0.5164
jenis kontrak: Two year dibanding kategori acuan ContractTwo year -1.3571 0.2574
tagihan tanpa kertas: Yes dibanding kategori acuan PaperlessBillingYes 0.3424 1.4083
metode pembayaran: Credit card (automatic) dibanding kategori acuan PaymentMethodCredit card (automatic) -0.0878 0.9160
metode pembayaran: Electronic check dibanding kategori acuan PaymentMethodElectronic check 0.3045 1.3559
metode pembayaran: Mailed check dibanding kategori acuan PaymentMethodMailed check -0.0576 0.9440
biaya bulanan pelanggan MonthlyCharges -0.0403 0.9605

6.7 Koefisien dan Odds Ratio Model Biner

tabel_biner <- broom::tidy(model_biner) %>%
  dplyr::mutate(
    Odds_Ratio       = exp(estimate),
    Variabel_Mudah   = purrr::map_chr(term, label_biner),
    Keputusan_5persen= ifelse(p.value < 0.05, "Signifikan", "Tidak signifikan")
  ) %>%
  dplyr::arrange(p.value)

knitr::kable(
  tabel_biner %>%
    dplyr::select(Variabel_Mudah, term, estimate, std.error,
                  statistic, p.value, Odds_Ratio, Keputusan_5persen) %>%
    head(25),
  digits  = 4,
  caption = "Koefisien dan Odds Ratio Regresi Logistik Biner"
)
Koefisien dan Odds Ratio Regresi Logistik Biner
Variabel_Mudah term estimate std.error statistic p.value Odds_Ratio Keputusan_5persen
durasi berlangganan (bulan) tenure -0.0606 0.0062 -9.7163 0.0000 0.9412 Signifikan
jenis kontrak: Two year dibanding kategori acuan ContractTwo year -1.3571 0.1764 -7.6914 0.0000 0.2574 Signifikan
jenis kontrak: One year dibanding kategori acuan ContractOne year -0.6608 0.1076 -6.1421 0.0000 0.5164 Signifikan
total biaya kumulatif pelanggan TotalCharges 0.0003 0.0001 4.6574 0.0000 1.0003 Signifikan
tagihan tanpa kertas: Yes dibanding kategori acuan PaperlessBillingYes 0.3424 0.0745 4.5956 0.0000 1.4083 Signifikan
metode pembayaran: Electronic check dibanding kategori acuan PaymentMethodElectronic check 0.3045 0.0945 3.2220 0.0013 1.3559 Signifikan
status lansia (1 = Ya, 0 = Tidak) SeniorCitizen 0.2168 0.0845 2.5645 0.0103 1.2421 Signifikan
multiple lines: Yes dibanding kategori acuan MultipleLinesYes 0.4484 0.1773 2.5297 0.0114 1.5658 Signifikan
layanan internet: No dibanding kategori acuan InternetServiceNo -1.7863 0.8073 -2.2128 0.0269 0.1676 Signifikan
layanan internet: Fiber optic dibanding kategori acuan InternetServiceFiber optic 1.7475 0.7981 2.1896 0.0286 5.7401 Signifikan
layanan streaming film: Yes dibanding kategori acuan StreamingMoviesYes 0.5993 0.3267 1.8345 0.0666 1.8208 Tidak signifikan
layanan streaming TV: Yes dibanding kategori acuan StreamingTVYes 0.5905 0.3263 1.8097 0.0703 1.8049 Tidak signifikan
status tanggungan: Yes dibanding kategori acuan DependentsYes -0.1485 0.0897 -1.6548 0.0980 0.8620 Tidak signifikan
(Intercept) (Intercept) 1.1653 0.8151 1.4296 0.1528 3.2068 Tidak signifikan
biaya bulanan pelanggan MonthlyCharges -0.0403 0.0318 -1.2705 0.2039 0.9605 Tidak signifikan
keamanan online: Yes dibanding kategori acuan OnlineSecurityYes -0.2054 0.1787 -1.1496 0.2503 0.8143 Tidak signifikan
dukungan teknis: Yes dibanding kategori acuan TechSupportYes -0.1805 0.1806 -0.9994 0.3176 0.8349 Tidak signifikan
perlindungan perangkat: Yes dibanding kategori acuan DeviceProtectionYes 0.1474 0.1764 0.8356 0.4034 1.1588 Tidak signifikan
metode pembayaran: Credit card (automatic) dibanding kategori acuan PaymentMethodCredit card (automatic) -0.0878 0.1141 -0.7696 0.4416 0.9160 Tidak signifikan
metode pembayaran: Mailed check dibanding kategori acuan PaymentMethodMailed check -0.0576 0.1149 -0.5011 0.6163 0.9440 Tidak signifikan
jenis kelamin: Male dibanding kategori acuan genderMale -0.0218 0.0648 -0.3369 0.7362 0.9784 Tidak signifikan
layanan telepon: Yes dibanding kategori acuan PhoneServiceYes 0.1715 0.6487 0.2643 0.7915 1.1870 Tidak signifikan
pencadangan online: Yes dibanding kategori acuan OnlineBackupYes 0.0260 0.1754 0.1485 0.8820 1.0264 Tidak signifikan
status pasangan: Yes dibanding kategori acuan PartnerYes -0.0004 0.0778 -0.0049 0.9961 0.9996 Tidak signifikan
multiple lines: No phone service dibanding kategori acuan MultipleLinesNo phone service NA NA NA NA NA NA

6.8 Evaluasi Model Biner

data_biner$prob_churn <- predict(model_biner, type = "response")
data_biner$pred_churn <- ifelse(data_biner$prob_churn >= 0.5, 1, 0)

tabel_klasifikasi_biner <- table(
  Prediksi = factor(data_biner$pred_churn, levels = c(0, 1),
                    labels = c("Tidak churn", "Churn")),
  Aktual   = data_biner$Churn_label
)
tabel_klasifikasi_biner
##              Aktual
## Prediksi        No  Yes
##   Tidak churn 4639  833
##   Churn        524 1036
akurasi_biner <- mean(data_biner$pred_churn == data_biner$Churn_num)
akurasi_biner
## [1] 0.807025
AIC_biner <- AIC(model_biner)
AIC_biner
## [1] 5874.272

6.9 Visualisasi Probabilitas Prediksi Model Biner

Visualisasi ini digunakan untuk melihat sebaran probabilitas churn yang dihasilkan oleh model. Probabilitas yang lebih tinggi menunjukkan pelanggan yang diprediksi memiliki risiko churn lebih besar.

ggplot(data_biner, aes(x = prob_churn)) +
  geom_histogram(bins = 40, fill = "#2f7f73", color = "white") +
  labs(
    title = "Distribusi Probabilitas Prediksi Churn",
    x     = "Probabilitas prediksi churn",
    y     = "Jumlah pelanggan"
  ) +
  theme_minimal()

ggplot(data_biner, aes(x = Churn_label, y = prob_churn, fill = Churn_label)) +
  geom_boxplot() +
  scale_fill_manual(values = c("No" = "#2f7f73", "Yes" = "#b84f5a")) +
  labs(
    title = "Probabilitas Prediksi Churn menurut Status Aktual",
    x     = "Status aktual",
    y     = "Probabilitas prediksi churn",
    fill  = "Status churn"
  ) +
  theme_minimal()

Grafik pertama menunjukkan sebaran probabilitas churn hasil model. Grafik kedua membandingkan probabilitas prediksi antara pelanggan yang benar-benar churn dan tidak churn. Jika model cukup membedakan risiko, maka kelompok yang aktual churn seharusnya memiliki probabilitas prediksi yang lebih tinggi.

6.10 Evaluasi Tambahan: ROC Curve dan AUC

ROC Curve digunakan untuk melihat kemampuan model membedakan pelanggan churn dan tidak churn pada berbagai threshold. AUC mendekati 1 menunjukkan kemampuan pemisahan yang lebih baik, sedangkan AUC mendekati 0,5 menunjukkan kemampuan model mendekati tebakan acak.

roc_biner <- pROC::roc(
  response  = data_biner$Churn_num,
  predictor = data_biner$prob_churn,
  quiet     = TRUE
)

auc_biner <- pROC::auc(roc_biner)
auc_biner
## Area under the curve: 0.8482
plot(
  roc_biner,
  main = paste("ROC Curve Model Biner | AUC =", round(as.numeric(auc_biner), 3))
)

AUC menunjukkan kemampuan model membedakan pelanggan yang churn dan tidak churn. Semakin besar AUC, semakin baik kemampuan diskriminasi model. Namun AUC tidak menjelaskan arah pengaruh variabel; arah pengaruh tetap dibaca melalui koefisien dan odds ratio.

6.11 Catatan Asumsi Regresi Logistik Biner

Regresi logistik biner digunakan karena variabel respon hanya memiliki dua kategori, yaitu churn dan tidak churn. Model ini mengasumsikan bahwa observasi bersifat independen, tidak terdapat multikolinearitas berat antar prediktor, tidak terjadi perfect separation, dan prediktor numerik memiliki hubungan yang cukup linear terhadap log-odds churn.

Dalam analisis churn, jumlah pelanggan yang churn sering kali lebih kecil daripada pelanggan yang tidak churn (kelas tidak seimbang). Oleh karena itu, akurasi tidak boleh menjadi satu-satunya ukuran evaluasi. Confusion matrix, probabilitas prediksi, dan AUC perlu dibaca bersama dengan interpretasi odds ratio.

6.12 Interpretasi dan Kesimpulan Model Biner

sig_biner <- tabel_biner %>%
  dplyr::filter(term != "(Intercept)", p.value < 0.05) %>%
  dplyr::arrange(p.value)

cat("Pada taraf signifikansi 5%, model biner menemukan ", nrow(sig_biner),
    " prediktor yang berhubungan signifikan dengan peluang churn pelanggan.\n\n", sep = "")

Pada taraf signifikansi 5%, model biner menemukan 10 prediktor yang berhubungan signifikan dengan peluang churn pelanggan.

if (nrow(sig_biner) > 0) {
  cat("Prediktor paling penting berdasarkan p-value terkecil adalah:\n\n")
  top_biner <- head(sig_biner, 8)
  for (i in seq_len(nrow(top_biner))) {
    row  <- top_biner[i, ]
    arah <- ifelse(row$Odds_Ratio > 1,
                   "meningkatkan odds churn pelanggan",
                   "menurunkan odds churn pelanggan")
    cat("- **", row$Variabel_Mudah, "** memiliki OR = ", fmt_num(row$Odds_Ratio, 3),
        " dan p-value = ", fmt_p(row$p.value), ", sehingga variabel ini ", arah,
        ", dengan asumsi variabel lain tetap.\n", sep = "")
  }
  cat("\n")
}

Prediktor paling penting berdasarkan p-value terkecil adalah:

  • durasi berlangganan (bulan) memiliki OR = 0.941 dan p-value = < 0.001, sehingga variabel ini menurunkan odds churn pelanggan, dengan asumsi variabel lain tetap.
  • jenis kontrak: Two year dibanding kategori acuan memiliki OR = 0.257 dan p-value = < 0.001, sehingga variabel ini menurunkan odds churn pelanggan, dengan asumsi variabel lain tetap.
  • jenis kontrak: One year dibanding kategori acuan memiliki OR = 0.516 dan p-value = < 0.001, sehingga variabel ini menurunkan odds churn pelanggan, dengan asumsi variabel lain tetap.
  • total biaya kumulatif pelanggan memiliki OR = 1.000 dan p-value = < 0.001, sehingga variabel ini meningkatkan odds churn pelanggan, dengan asumsi variabel lain tetap.
  • tagihan tanpa kertas: Yes dibanding kategori acuan memiliki OR = 1.408 dan p-value = < 0.001, sehingga variabel ini meningkatkan odds churn pelanggan, dengan asumsi variabel lain tetap.
  • metode pembayaran: Electronic check dibanding kategori acuan memiliki OR = 1.356 dan p-value = 0.001, sehingga variabel ini meningkatkan odds churn pelanggan, dengan asumsi variabel lain tetap.
  • status lansia (1 = Ya, 0 = Tidak) memiliki OR = 1.242 dan p-value = 0.010, sehingga variabel ini meningkatkan odds churn pelanggan, dengan asumsi variabel lain tetap.
  • multiple lines: Yes dibanding kategori acuan memiliki OR = 1.566 dan p-value = 0.011, sehingga variabel ini meningkatkan odds churn pelanggan, dengan asumsi variabel lain tetap.
cat("Akurasi klasifikasi model biner adalah **", scales::percent(akurasi_biner),
    "**. AIC model adalah **", fmt_num(AIC_biner, 2), "**.\n\n", sep = "")

Akurasi klasifikasi model biner adalah 81%. AIC model adalah 5,874.27.

cat("**Kesimpulan substantif:** jenis kontrak, lama berlangganan, dan jenis layanan internet merupakan faktor utama yang menjelaskan peluang churn pelanggan telekomunikasi. Pelanggan dengan kontrak bulanan dan berlangganan lebih singkat cenderung memiliki risiko churn yang lebih tinggi.\n")

Kesimpulan substantif: jenis kontrak, lama berlangganan, dan jenis layanan internet merupakan faktor utama yang menjelaskan peluang churn pelanggan telekomunikasi. Pelanggan dengan kontrak bulanan dan berlangganan lebih singkat cenderung memiliki risiko churn yang lebih tinggi.

7 Bagian II — Regresi Logistik Ordinal

7.1 Tujuan Analisis

Tujuan analisis ini adalah mengetahui faktor-faktor yang berhubungan dengan tingkat performa akademik siswa yang diukur melalui nilai akhir (G3) pada mata pelajaran bahasa Portugis. Nilai G3 dikelompokkan menjadi tiga kategori ordinal: Rendah, Sedang, dan Tinggi.

7.2 Variabel yang Digunakan

Variabel respon adalah Performa_G3 dengan urutan:

\[ \text{Rendah} < \text{Sedang} < \text{Tinggi} \]

Variabel yang Dipakai pada Model Ordinal
Variabel Peran Alasan_dipakai
Performa_G3 Respon Nilai akhir G3 dikelompokkan menjadi Rendah/Sedang/Tinggi
studytime Prediktor Waktu belajar berkaitan langsung dengan capaian nilai
failures Prediktor Riwayat kegagalan akademik sebelumnya
Dalc Prediktor Konsumsi alkohol hari kerja (skala 1-5)
Walc Prediktor Konsumsi alkohol akhir pekan (skala 1-5)
absences Prediktor Banyaknya ketidakhadiran dapat menurunkan capaian
famrel Prediktor Hubungan keluarga yang baik mendukung belajar
goout Prediktor Frekuensi bermain bersama teman
health Prediktor Status kesehatan siswa
Medu Prediktor Pendidikan ibu berkaitan dengan dukungan belajar di rumah
Fedu Prediktor Pendidikan ayah berkaitan dengan dukungan belajar di rumah
age Prediktor Usia siswa
internet Prediktor Akses internet di rumah
higher Prediktor Motivasi melanjutkan pendidikan tinggi
traveltime Prediktor Waktu perjalanan ke sekolah

Variabel G1 dan G2 tidak digunakan dalam model karena merupakan nilai antara yang sangat berkorelasi dengan G3. Penggunaan keduanya akan menyebabkan data leakage karena nilai akhir G3 sangat ditentukan oleh G1 dan G2. Selain itu, variabel school, sex, dan variabel kategorik lainnya digunakan sebagai prediktor tambahan.

7.3 Cleaning Data Ordinal

data_ordinal <- data_ordinal_raw %>%
  dplyr::mutate(
    Performa_G3 = dplyr::case_when(
      G3 <= 9  ~ "Rendah",
      G3 <= 14 ~ "Sedang",
      G3 >= 15 ~ "Tinggi"
    ),
    Performa_G3 = ordered(Performa_G3, levels = c("Rendah", "Sedang", "Tinggi")),
    school      = factor(school),
    sex         = factor(sex),
    address     = factor(address),
    famsize     = factor(famsize),
    Pstatus     = factor(Pstatus),
    Mjob        = factor(Mjob),
    Fjob        = factor(Fjob),
    reason      = factor(reason),
    guardian    = factor(guardian),
    schoolsup   = factor(schoolsup),
    famsup      = factor(famsup),
    paid        = factor(paid),
    activities  = factor(activities),
    nursery     = factor(nursery),
    higher      = factor(higher),
    internet    = factor(internet),
    romantic    = factor(romantic)
  ) %>%
  tidyr::drop_na()

knitr::kable(
  head(data_ordinal %>%
         dplyr::select(school, sex, age, studytime, failures, Dalc, Walc,
                       absences, G3, Performa_G3), 6),
  caption = "Cuplikan Data Ordinal Setelah Cleaning"
)
Cuplikan Data Ordinal Setelah Cleaning
school sex age studytime failures Dalc Walc absences G3 Performa_G3
GP F 18 2 0 1 1 4 11 Sedang
GP F 17 2 0 1 1 2 11 Sedang
GP F 15 2 0 2 3 6 12 Sedang
GP F 15 3 0 1 1 0 14 Sedang
GP F 16 2 0 1 2 0 13 Sedang
GP M 16 2 0 1 2 6 13 Sedang

Nilai G3 dikelompokkan menjadi tiga kategori ordinal: Rendah (0–9), Sedang (10–14), dan Tinggi (15–20). Pengelompokan ini mengacu pada konvensi sistem penilaian Portugis dan bertujuan agar model ordinal dapat diaplikasikan dengan tiga tingkatan yang bermakna secara substantif.

7.4 Analisis Deskriptif Data Ordinal

dist_ordinal <- data_ordinal %>%
  dplyr::count(Performa_G3) %>%
  dplyr::mutate(Persen = scales::percent(n / sum(n)))
knitr::kable(dist_ordinal, caption = "Distribusi Performa G3 Siswa")
Distribusi Performa G3 Siswa
Performa_G3 n Persen
Rendah 100 15.4%
Sedang 418 64.4%
Tinggi 131 20.2%
ringkasan_ordinal <- data_ordinal %>%
  dplyr::summarise(
    Rata_rata_G3       = mean(G3),
    Median_G3          = median(G3),
    Rata_rata_studytime= mean(studytime),
    Rata_rata_absences = mean(absences),
    Rata_rata_failures = mean(failures)
  )
knitr::kable(ringkasan_ordinal, caption = "Ringkasan Prediktor Utama Data Ordinal")
Ringkasan Prediktor Utama Data Ordinal
Rata_rata_G3 Median_G3 Rata_rata_studytime Rata_rata_absences Rata_rata_failures
11.90601 12 1.930663 3.659476 0.2218798
ggplot(data_ordinal, aes(x = Performa_G3, fill = Performa_G3)) +
  geom_bar() +
  scale_fill_manual(values = c("Rendah" = "#b84f5a", "Sedang" = "#d18b2f", "Tinggi" = "#2f7f73")) +
  labs(
    title = "Distribusi Performa Akademik Siswa (G3)",
    x     = "Kategori performa",
    y     = "Jumlah siswa",
    fill  = "Performa"
  ) +
  theme_minimal()

7.5 Pembentukan Model Regresi Logistik Ordinal

Model ordinal dibentuk untuk menjelaskan faktor-faktor yang berhubungan dengan tingkat performa akademik siswa. Model ini menggunakan spesifikasi MASS::polr() sebagai berikut:

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

Dalam penelitian ini, variabel respon adalah Performa_G3, dengan urutan:

\[ \text{Rendah} < \text{Sedang} < \text{Tinggi} \]

Karena model menggunakan bentuk \(\alpha_j - x\beta\), maka koefisien positif menunjukkan penurunan peluang kumulatif berada pada kategori rendah. Dengan kata lain, koefisien positif menunjukkan peningkatan kecenderungan siswa untuk bergeser ke tingkat performa yang lebih tinggi.

model_ordinal <- MASS::polr(
  Performa_G3 ~ school + sex + age + address + famsize + Pstatus +
    Medu + Fedu + Mjob + Fjob + reason + guardian +
    traveltime + studytime + failures + schoolsup + famsup + paid +
    activities + nursery + higher + internet + romantic +
    famrel + freetime + goout + Dalc + Walc + health + absences,
  data = data_ordinal,
  Hess = TRUE
)

summary(model_ordinal)
## Call:
## MASS::polr(formula = Performa_G3 ~ school + sex + age + address + 
##     famsize + Pstatus + Medu + Fedu + Mjob + Fjob + reason + 
##     guardian + traveltime + studytime + failures + schoolsup + 
##     famsup + paid + activities + nursery + higher + internet + 
##     romantic + famrel + freetime + goout + Dalc + Walc + health + 
##     absences, data = data_ordinal, Hess = TRUE)
## 
## Coefficients:
##                      Value Std. Error  t value
## schoolMS         -1.135536    0.23216 -4.89120
## sexM             -0.446014    0.20402 -2.18617
## age               0.280442    0.08380  3.34650
## addressU          0.009603    0.21598  0.04447
## famsizeLE3        0.295601    0.20013  1.47705
## PstatusT         -0.014614    0.28619 -0.05107
## Medu              0.172790    0.12616  1.36960
## Fedu              0.066272    0.11440  0.57928
## Mjobhealth        0.441884    0.43667  1.01195
## Mjobother         0.026157    0.24728  0.10578
## Mjobservices      0.247434    0.30451  0.81257
## Mjobteacher       0.218389    0.40863  0.53444
## Fjobhealth       -0.603004    0.62190 -0.96962
## Fjobother        -0.328684    0.37268 -0.88194
## Fjobservices     -0.595829    0.39149 -1.52194
## Fjobteacher       0.134257    0.54458  0.24653
## reasonhome        0.087483    0.23226  0.37666
## reasonother      -0.097900    0.30474 -0.32126
## reasonreputation  0.318822    0.24219  1.31639
## guardianmother   -0.342715    0.21594 -1.58709
## guardianother    -0.228469    0.43264 -0.52809
## traveltime        0.127330    0.13257  0.96044
## studytime         0.255165    0.11373  2.24355
## failures         -1.024736    0.18080 -5.66772
## schoolsupyes     -1.076297    0.30649 -3.51163
## famsupyes        -0.087918    0.18617 -0.47223
## paidyes          -0.619844    0.38686 -1.60223
## activitiesyes     0.087654    0.18229  0.48085
## nurseryyes       -0.189096    0.22327 -0.84693
## higheryes         1.492243    0.31799  4.69273
## internetyes       0.099099    0.22815  0.43437
## romanticyes      -0.119833    0.18764 -0.63865
## famrel            0.098758    0.09529  1.03643
## freetime         -0.026837    0.09290 -0.28889
## goout            -0.055981    0.08845 -0.63289
## Dalc             -0.142805    0.12732 -1.12161
## Walc             -0.039912    0.09733 -0.41009
## health           -0.110034    0.06371 -1.72701
## absences         -0.059224    0.02117 -2.79714
## 
## Intercepts:
##               Value   Std. Error t value
## Rendah|Sedang  2.7716  1.6248     1.7059
## Sedang|Tinggi  6.9768  1.6527     4.2214
## 
## Residual Deviance: 927.9076 
## AIC: 1009.908

7.6 Persamaan Model Hasil Estimasi Ordinal

Chunk model-ordinal di atas sudah membangun model regresi logistik ordinal menggunakan MASS::polr(). Bentuk hasil estimasinya adalah:

\[ \log\left(\frac{\widehat{P}(\text{Performa}_i \leq j)}{\widehat{P}(\text{Performa}_i > j)}\right) = \hat{\alpha}_j - \widehat{\eta_i} \]

Dengan:

\[ \widehat{\eta_i} = \hat{\beta}_1X_{1i} + \hat{\beta}_2X_{2i} + \cdots + \hat{\beta}_pX_{pi} \]

Karena MASS::polr() memakai bentuk alpha_j - x beta, maka exp(beta) dibaca sebagai odds menuju performa yang lebih tinggi.

Parameter yang Membentuk Persamaan Regresi Logistik Ordinal
Parameter Estimate Std_Error t_value Jenis Odds_Ratio_Higher_Perf
schoolMS -1.1355 0.2322 -4.8912 Koefisien prediktor 0.3212
sexM -0.4460 0.2040 -2.1862 Koefisien prediktor 0.6402
age 0.2804 0.0838 3.3465 Koefisien prediktor 1.3237
addressU 0.0096 0.2160 0.0445 Koefisien prediktor 1.0096
famsizeLE3 0.2956 0.2001 1.4770 Koefisien prediktor 1.3439
PstatusT -0.0146 0.2862 -0.0511 Koefisien prediktor 0.9855
Medu 0.1728 0.1262 1.3696 Koefisien prediktor 1.1886
Fedu 0.0663 0.1144 0.5793 Koefisien prediktor 1.0685
Mjobhealth 0.4419 0.4367 1.0120 Koefisien prediktor 1.5556
Mjobother 0.0262 0.2473 0.1058 Koefisien prediktor 1.0265
Mjobservices 0.2474 0.3045 0.8126 Koefisien prediktor 1.2807
Mjobteacher 0.2184 0.4086 0.5344 Koefisien prediktor 1.2441
Fjobhealth -0.6030 0.6219 -0.9696 Koefisien prediktor 0.5472
Fjobother -0.3287 0.3727 -0.8819 Koefisien prediktor 0.7199
Fjobservices -0.5958 0.3915 -1.5219 Koefisien prediktor 0.5511
Fjobteacher 0.1343 0.5446 0.2465 Koefisien prediktor 1.1437
reasonhome 0.0875 0.2323 0.3767 Koefisien prediktor 1.0914
reasonother -0.0979 0.3047 -0.3213 Koefisien prediktor 0.9067
reasonreputation 0.3188 0.2422 1.3164 Koefisien prediktor 1.3755
guardianmother -0.3427 0.2159 -1.5871 Koefisien prediktor 0.7098
guardianother -0.2285 0.4326 -0.5281 Koefisien prediktor 0.7958
traveltime 0.1273 0.1326 0.9604 Koefisien prediktor 1.1358
studytime 0.2552 0.1137 2.2435 Koefisien prediktor 1.2907
failures -1.0247 0.1808 -5.6677 Koefisien prediktor 0.3589
schoolsupyes -1.0763 0.3065 -3.5116 Koefisien prediktor 0.3409
famsupyes -0.0879 0.1862 -0.4722 Koefisien prediktor 0.9158
paidyes -0.6198 0.3869 -1.6022 Koefisien prediktor 0.5380
activitiesyes 0.0877 0.1823 0.4808 Koefisien prediktor 1.0916
nurseryyes -0.1891 0.2233 -0.8469 Koefisien prediktor 0.8277
higheryes 1.4922 0.3180 4.6927 Koefisien prediktor 4.4471
internetyes 0.0991 0.2281 0.4344 Koefisien prediktor 1.1042
romanticyes -0.1198 0.1876 -0.6386 Koefisien prediktor 0.8871
famrel 0.0988 0.0953 1.0364 Koefisien prediktor 1.1038
freetime -0.0268 0.0929 -0.2889 Koefisien prediktor 0.9735
goout -0.0560 0.0885 -0.6329 Koefisien prediktor 0.9456
Dalc -0.1428 0.1273 -1.1216 Koefisien prediktor 0.8669
Walc -0.0399 0.0973 -0.4101 Koefisien prediktor 0.9609
health -0.1100 0.0637 -1.7270 Koefisien prediktor 0.8958
absences -0.0592 0.0212 -2.7971 Koefisien prediktor 0.9425
Rendah|Sedang 2.7716 1.6248 1.7059 Cutpoint NA
Sedang|Tinggi 6.9768 1.6527 4.2214 Cutpoint NA

7.7 Koefisien dan Odds Ratio Model Ordinal

Pada tabel berikut, Odds_Ratio_Higher_Perf = exp(estimate) dibaca sebagai odds untuk berada pada kategori performa yang lebih tinggi. Ini sesuai dengan parameterisasi MASS::polr() yang berbentuk:

\[ \operatorname{logit}\{P(Y \leq j)\} = \alpha_j - x\beta \]

coef_ordinal <- coef(summary(model_ordinal)) %>%
  as.data.frame()
coef_ordinal$term <- rownames(coef_ordinal)
names(coef_ordinal)[1:3] <- c("estimate", "std_error", "t_value")

coef_ordinal <- coef_ordinal %>%
  dplyr::mutate(
    p_value                = 2 * pnorm(abs(t_value), lower.tail = FALSE),
    Odds_Ratio_Higher_Perf = exp(estimate),
    Odds_Ratio_Lower_Cumul = exp(-estimate),
    Jenis                  = ifelse(stringr::str_detect(term, "\\|"), "Cutpoint", "Prediktor"),
    Variabel_Mudah         = purrr::map_chr(term, label_ordinal),
    Keputusan_5persen      = ifelse(p_value < 0.05, "Signifikan", "Tidak signifikan")
  ) %>%
  dplyr::relocate(term, Jenis)

knitr::kable(
  coef_ordinal %>%
    dplyr::filter(Jenis == "Prediktor") %>%
    dplyr::arrange(p_value) %>%
    dplyr::select(Variabel_Mudah, term, estimate, std_error, t_value,
                  p_value, Odds_Ratio_Higher_Perf, Odds_Ratio_Lower_Cumul,
                  Keputusan_5persen) %>%
    head(25),
  digits  = 4,
  caption = "Koefisien dan Odds Ratio Regresi Ordinal MASS::polr"
)
Koefisien dan Odds Ratio Regresi Ordinal MASS::polr
Variabel_Mudah term estimate std_error t_value p_value Odds_Ratio_Higher_Perf Odds_Ratio_Lower_Cumul Keputusan_5persen
failures jumlah kegagalan kelas sebelumnya failures -1.0247 0.1808 -5.6677 0.0000 0.3589 2.7864 Signifikan
schoolMS sekolah: MS dibanding kategori acuan schoolMS -1.1355 0.2322 -4.8912 0.0000 0.3212 3.1128 Signifikan
higheryes ingin melanjutkan pendidikan tinggi: yes dibanding acuan higheryes 1.4922 0.3180 4.6927 0.0000 4.4471 0.2249 Signifikan
schoolsupyes sekolah: supyes dibanding kategori acuan schoolsupyes -1.0763 0.3065 -3.5116 0.0004 0.3409 2.9338 Signifikan
age usia siswa age 0.2804 0.0838 3.3465 0.0008 1.3237 0.7554 Signifikan
absences jumlah absensi/ketidakhadiran sekolah absences -0.0592 0.0212 -2.7971 0.0052 0.9425 1.0610 Signifikan
studytime waktu belajar mingguan (skala 1-4) studytime 0.2552 0.1137 2.2435 0.0249 1.2907 0.7748 Signifikan
sexM jenis kelamin siswa: M dibanding kategori acuan sexM -0.4460 0.2040 -2.1862 0.0288 0.6402 1.5621 Signifikan
health status kesehatan saat ini health -0.1100 0.0637 -1.7270 0.0842 0.8958 1.1163 Tidak signifikan
paidyes les berbayar: yes dibanding acuan paidyes -0.6198 0.3869 -1.6022 0.1091 0.5380 1.8586 Tidak signifikan
guardianmother wali siswa: mother dibanding acuan guardianmother -0.3427 0.2159 -1.5871 0.1125 0.7098 1.4088 Tidak signifikan
Fjobservices pekerjaan ayah: services dibanding kategori acuan Fjobservices -0.5958 0.3915 -1.5219 0.1280 0.5511 1.8145 Tidak signifikan
famsizeLE3 ukuran keluarga: LE3 dibanding acuan famsizeLE3 0.2956 0.2001 1.4770 0.1397 1.3439 0.7441 Tidak signifikan
Medu tingkat pendidikan ibu Medu 0.1728 0.1262 1.3696 0.1708 1.1886 0.8413 Tidak signifikan
reasonreputation alasan pilih sekolah: reputation dibanding acuan reasonreputation 0.3188 0.2422 1.3164 0.1880 1.3755 0.7270 Tidak signifikan
Dalc konsumsi alkohol hari kerja (skala 1-5) Dalc -0.1428 0.1273 -1.1216 0.2620 0.8669 1.1535 Tidak signifikan
famrel kualitas hubungan keluarga famrel 0.0988 0.0953 1.0364 0.3000 1.1038 0.9060 Tidak signifikan
Mjobhealth pekerjaan ibu: health dibanding kategori acuan Mjobhealth 0.4419 0.4367 1.0120 0.3116 1.5556 0.6428 Tidak signifikan
Fjobhealth pekerjaan ayah: health dibanding kategori acuan Fjobhealth -0.6030 0.6219 -0.9696 0.3322 0.5472 1.8276 Tidak signifikan
traveltime waktu perjalanan ke sekolah traveltime 0.1273 0.1326 0.9604 0.3368 1.1358 0.8804 Tidak signifikan
Fjobother pekerjaan ayah: other dibanding kategori acuan Fjobother -0.3287 0.3727 -0.8819 0.3778 0.7199 1.3891 Tidak signifikan
nurseryyes pernah ikut taman kanak-kanak: yes dibanding acuan nurseryyes -0.1891 0.2233 -0.8469 0.3970 0.8277 1.2082 Tidak signifikan
Mjobservices pekerjaan ibu: services dibanding kategori acuan Mjobservices 0.2474 0.3045 0.8126 0.4165 1.2807 0.7808 Tidak signifikan
romanticyes memiliki hubungan romantis: yes dibanding acuan romanticyes -0.1198 0.1876 -0.6386 0.5231 0.8871 1.1273 Tidak signifikan
goout frekuensi pergi bermain bersama teman goout -0.0560 0.0885 -0.6329 0.5268 0.9456 1.0576 Tidak signifikan

Pada model ordinal dengan MASS::polr(), nilai Odds_Ratio_Higher_Perf = exp(beta) digunakan untuk membaca kecenderungan menuju performa yang lebih tinggi. Karena performa diurutkan sebagai Rendah < Sedang < Tinggi, maka OR di atas 1 berarti variabel tersebut berkaitan dengan performa akademik yang lebih baik, sedangkan OR di bawah 1 berarti berkaitan dengan performa yang lebih rendah.

7.8 Evaluasi Model Ordinal

pred_ordinal     <- predict(model_ordinal, newdata = data_ordinal, type = "class")
perf_levels      <- levels(data_ordinal$Performa_G3)
pred_ordinal_f   <- factor(as.character(pred_ordinal), levels = perf_levels, ordered = TRUE)
aktual_ordinal_f <- factor(as.character(data_ordinal$Performa_G3), levels = perf_levels, ordered = TRUE)

conf_ordinal <- table(Prediksi = pred_ordinal_f, Aktual = aktual_ordinal_f)
conf_ordinal
##         Aktual
## Prediksi Rendah Sedang Tinggi
##   Rendah     30     19      0
##   Sedang     69    383    105
##   Tinggi      1     16     26
akurasi_ordinal <- mean(as.character(pred_ordinal_f) == as.character(aktual_ordinal_f), na.rm = TRUE)
akurasi_ordinal
## [1] 0.6764253
AIC_ordinal <- AIC(model_ordinal)
AIC_ordinal
## [1] 1009.908
prob_ordinal <- predict(model_ordinal, type = "probs")
knitr::kable(
  head(prob_ordinal, 6),
  digits  = 4,
  caption = "Cuplikan Prediksi Probabilitas untuk Model Ordinal"
)
Cuplikan Prediksi Probabilitas untuk Model Ordinal
Rendah Sedang Tinggi
0.0334 0.6648 0.3018
0.0252 0.6087 0.3661
0.2191 0.7304 0.0505
0.0318 0.6556 0.3126
0.0294 0.6409 0.3296
0.0284 0.6340 0.3375

7.9 Visualisasi Probabilitas Prediksi Model Ordinal

Visualisasi ini digunakan untuk melihat bagaimana peluang masing-masing kategori performa berubah ketika absences (jumlah absensi) meningkat. Variabel lain ditahan pada nilai median atau kategori paling umum.

grid_ordinal <- data.frame(
  absences   = seq(
    quantile(data_ordinal$absences, 0.05, na.rm = TRUE),
    quantile(data_ordinal$absences, 0.95, na.rm = TRUE),
    length.out = 100
  ),
  school     = factor(modus(data_ordinal$school),  levels = levels(data_ordinal$school)),
  sex        = factor(modus(data_ordinal$sex),     levels = levels(data_ordinal$sex)),
  age        = median(data_ordinal$age, na.rm = TRUE),
  address    = factor(modus(data_ordinal$address), levels = levels(data_ordinal$address)),
  famsize    = factor(modus(data_ordinal$famsize), levels = levels(data_ordinal$famsize)),
  Pstatus    = factor(modus(data_ordinal$Pstatus), levels = levels(data_ordinal$Pstatus)),
  Medu       = median(data_ordinal$Medu,      na.rm = TRUE),
  Fedu       = median(data_ordinal$Fedu,      na.rm = TRUE),
  Mjob       = factor(modus(data_ordinal$Mjob),    levels = levels(data_ordinal$Mjob)),
  Fjob       = factor(modus(data_ordinal$Fjob),    levels = levels(data_ordinal$Fjob)),
  reason     = factor(modus(data_ordinal$reason),  levels = levels(data_ordinal$reason)),
  guardian   = factor(modus(data_ordinal$guardian),levels = levels(data_ordinal$guardian)),
  traveltime = median(data_ordinal$traveltime, na.rm = TRUE),
  studytime  = median(data_ordinal$studytime,  na.rm = TRUE),
  failures   = 0,
  schoolsup  = factor(modus(data_ordinal$schoolsup), levels = levels(data_ordinal$schoolsup)),
  famsup     = factor(modus(data_ordinal$famsup),    levels = levels(data_ordinal$famsup)),
  paid       = factor(modus(data_ordinal$paid),      levels = levels(data_ordinal$paid)),
  activities = factor(modus(data_ordinal$activities),levels = levels(data_ordinal$activities)),
  nursery    = factor(modus(data_ordinal$nursery),   levels = levels(data_ordinal$nursery)),
  higher     = factor(modus(data_ordinal$higher),    levels = levels(data_ordinal$higher)),
  internet   = factor(modus(data_ordinal$internet),  levels = levels(data_ordinal$internet)),
  romantic   = factor(modus(data_ordinal$romantic),  levels = levels(data_ordinal$romantic)),
  famrel     = median(data_ordinal$famrel,  na.rm = TRUE),
  freetime   = median(data_ordinal$freetime,na.rm = TRUE),
  goout      = median(data_ordinal$goout,   na.rm = TRUE),
  Dalc       = median(data_ordinal$Dalc,    na.rm = TRUE),
  Walc       = median(data_ordinal$Walc,    na.rm = TRUE),
  health     = median(data_ordinal$health,  na.rm = TRUE)
)

prob_grid_ordinal <- predict(model_ordinal, newdata = grid_ordinal, type = "probs")

plot_prob_ordinal <- grid_ordinal %>%
  dplyr::bind_cols(as.data.frame(prob_grid_ordinal)) %>%
  tidyr::pivot_longer(
    cols      = levels(data_ordinal$Performa_G3),
    names_to  = "Performa_G3",
    values_to = "Probabilitas"
  ) %>%
  dplyr::mutate(
    Performa_G3 = factor(Performa_G3, levels = levels(data_ordinal$Performa_G3))
  )

ggplot(
  plot_prob_ordinal,
  aes(x = absences, y = Probabilitas, color = Performa_G3)
) +
  geom_line(linewidth = 1.1) +
  scale_color_manual(values = c("Rendah" = "#b84f5a", "Sedang" = "#d18b2f", "Tinggi" = "#2f7f73")) +
  scale_y_continuous(labels = scales::percent_format()) +
  labs(
    title    = "Prediksi Probabilitas Performa G3 Berdasarkan Jumlah Absensi",
    subtitle = "Variabel lain ditahan pada median atau kategori paling umum",
    x        = "Jumlah absensi",
    y        = "Probabilitas prediksi",
    color    = "Performa"
  ) +
  theme_minimal()

Grafik ini membantu melihat perubahan peluang kategori performa ketika jumlah absensi meningkat. Jika probabilitas kategori Rendah meningkat ketika absensi naik, maka model menunjukkan bahwa ketidakhadiran yang lebih banyak berkaitan dengan performa akademik yang lebih rendah.

7.10 Pemeriksaan Asumsi Proportional Odds

Asumsi utama regresi logistik ordinal adalah proportional odds atau parallel lines. Artinya, pengaruh prediktor diasumsikan sama pada setiap batas kumulatif performa.

cek_proportional_odds <- tryCatch({
  model_ordinal_clm <- ordinal::clm(
    Performa_G3 ~ school + sex + age + address + famsize + Pstatus +
      Medu + Fedu + Mjob + Fjob + reason + guardian +
      traveltime + studytime + failures + schoolsup + famsup + paid +
      activities + nursery + higher + internet + romantic +
      famrel + freetime + goout + Dalc + Walc + health + absences,
    data = data_ordinal,
    link = "logit"
  )

  print(summary(model_ordinal_clm))

  ordinal::nominal_test(model_ordinal_clm)
}, error = function(e) {
  message("Pemeriksaan proportional odds dengan nominal_test() tidak berhasil dijalankan: ", e$message)
  NULL
})
## formula: 
## Performa_G3 ~ school + sex + age + address + famsize + Pstatus + Medu + Fedu + Mjob + Fjob + reason + guardian + traveltime + studytime + failures + schoolsup + famsup + paid + activities + nursery + higher + internet + romantic + famrel + freetime + goout + Dalc + Walc + health + absences
## data:    data_ordinal
## 
##  link  threshold nobs logLik  AIC     niter max.grad cond.H 
##  logit flexible  649  -463.95 1009.91 6(0)  3.96e-07 2.8e+05
## 
## Coefficients:
##                   Estimate Std. Error z value     Pr(>|z|)    
## schoolMS         -1.135542   0.232158  -4.891 0.0000010020 ***
## sexM             -0.446020   0.204016  -2.186     0.028801 *  
## age               0.280465   0.083800   3.347     0.000817 ***
## addressU          0.009618   0.215975   0.045     0.964478    
## famsizeLE3        0.295612   0.200129   1.477     0.139647    
## PstatusT         -0.014650   0.286185  -0.051     0.959174    
## Medu              0.172778   0.126160   1.370     0.170840    
## Fedu              0.066286   0.114404   0.579     0.562314    
## Mjobhealth        0.441937   0.436664   1.012     0.311503    
## Mjobother         0.026205   0.247282   0.106     0.915605    
## Mjobservices      0.247434   0.304506   0.813     0.416462    
## Mjobteacher       0.218447   0.408628   0.535     0.592936    
## Fjobhealth       -0.603218   0.621898  -0.970     0.332065    
## Fjobother        -0.328801   0.372680  -0.882     0.377635    
## Fjobservices     -0.595932   0.391492  -1.522     0.127958    
## Fjobteacher       0.134110   0.544577   0.246     0.805477    
## reasonhome        0.087477   0.232255   0.377     0.706441    
## reasonother      -0.097894   0.304740  -0.321     0.748029    
## reasonreputation  0.318791   0.242194   1.316     0.188087    
## guardianmother   -0.342690   0.215939  -1.587     0.112518    
## guardianother    -0.228479   0.432635  -0.528     0.597423    
## traveltime        0.127357   0.132574   0.961     0.336729    
## studytime         0.255167   0.113733   2.244     0.024860 *  
## failures         -1.024713   0.180801  -5.668 0.0000000145 ***
## schoolsupyes     -1.076327   0.306494  -3.512     0.000445 ***
## famsupyes        -0.087886   0.186174  -0.472     0.636880    
## paidyes          -0.619860   0.386861  -1.602     0.109094    
## activitiesyes     0.087664   0.182290   0.481     0.630582    
## nurseryyes       -0.189084   0.223271  -0.847     0.397061    
## higheryes         1.492255   0.317990   4.693 0.0000026953 ***
## internetyes       0.099105   0.228145   0.434     0.664003    
## romanticyes      -0.119837   0.187635  -0.639     0.523036    
## famrel            0.098765   0.095286   1.037     0.299968    
## freetime         -0.026826   0.092897  -0.289     0.772752    
## goout            -0.055985   0.088454  -0.633     0.526778    
## Dalc             -0.142799   0.127321  -1.122     0.262046    
## Walc             -0.039907   0.097326  -0.410     0.681777    
## health           -0.110030   0.063713  -1.727     0.084176 .  
## absences         -0.059222   0.021173  -2.797     0.005158 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Threshold coefficients:
##               Estimate Std. Error z value
## Rendah|Sedang    2.772      1.625   1.706
## Sedang|Tinggi    6.977      1.653   4.222
cek_proportional_odds

Jika hasil nominal_test() menunjukkan p-value kecil pada prediktor tertentu, maka terdapat indikasi bahwa efek prediktor tersebut tidak sepenuhnya memenuhi asumsi proportional odds. Namun hasil ini tidak boleh dibaca secara mekanis; interpretasi tetap perlu mempertimbangkan tujuan analisis, teori pendidikan, ukuran sampel, dan kesederhanaan model.

Jika asumsi proportional odds tidak terpenuhi, maka model ordinal sederhana dapat menjadi terlalu kaku. Alternatif yang dapat dipertimbangkan adalah partial proportional odds model, generalized ordered logit model, atau regresi multinomial sebagai pembanding non-ordinal. Dalam laporan ini, model MASS::polr() tetap digunakan sebagai model utama karena tujuan analisis adalah memahami penerapan regresi logistik ordinal pada data performa akademik.

7.11 Interpretasi dan Kesimpulan Model Ordinal

sig_ordinal <- coef_ordinal %>%
  dplyr::filter(Jenis == "Prediktor", p_value < 0.05) %>%
  dplyr::arrange(p_value)

cat("Pada taraf signifikansi 5%, model ordinal menemukan ", nrow(sig_ordinal),
    " prediktor yang berhubungan signifikan dengan performa akademik siswa.\n\n", sep = "")

Pada taraf signifikansi 5%, model ordinal menemukan 8 prediktor yang berhubungan signifikan dengan performa akademik siswa.

if (nrow(sig_ordinal) > 0) {
  cat("Prediktor paling penting berdasarkan p-value terkecil adalah:\n\n")
  top_ordinal <- head(sig_ordinal, 8)
  for (i in seq_len(nrow(top_ordinal))) {
    row  <- top_ordinal[i, ]
    arah <- ifelse(row$Odds_Ratio_Higher_Perf > 1,
                   "meningkatkan kecenderungan siswa masuk ke performa yang lebih tinggi",
                   "menurunkan kecenderungan siswa masuk ke performa yang lebih tinggi")
    cat("- **", row$Variabel_Mudah, "** memiliki OR kategori lebih tinggi = ",
        fmt_num(row$Odds_Ratio_Higher_Perf, 3),
        " dan p-value = ", fmt_p(row$p_value), ", sehingga variabel ini ", arah,
        ", dengan asumsi variabel lain tetap. Dalam parameterisasi `polr`, ini berasal dari model `logit(P(Y <= j)) = alpha_j - x beta`.\n", sep = "")
  }
  cat("\n")
}

Prediktor paling penting berdasarkan p-value terkecil adalah:

  • jumlah kegagalan kelas sebelumnya memiliki OR kategori lebih tinggi = 0.359 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan siswa masuk ke performa yang lebih tinggi, dengan asumsi variabel lain tetap. Dalam parameterisasi polr, ini berasal dari model logit(P(Y <= j)) = alpha_j - x beta.
  • sekolah: MS dibanding kategori acuan memiliki OR kategori lebih tinggi = 0.321 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan siswa masuk ke performa yang lebih tinggi, dengan asumsi variabel lain tetap. Dalam parameterisasi polr, ini berasal dari model logit(P(Y <= j)) = alpha_j - x beta.
  • ingin melanjutkan pendidikan tinggi: yes dibanding acuan memiliki OR kategori lebih tinggi = 4.447 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan siswa masuk ke performa yang lebih tinggi, dengan asumsi variabel lain tetap. Dalam parameterisasi polr, ini berasal dari model logit(P(Y <= j)) = alpha_j - x beta.
  • sekolah: supyes dibanding kategori acuan memiliki OR kategori lebih tinggi = 0.341 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan siswa masuk ke performa yang lebih tinggi, dengan asumsi variabel lain tetap. Dalam parameterisasi polr, ini berasal dari model logit(P(Y <= j)) = alpha_j - x beta.
  • usia siswa memiliki OR kategori lebih tinggi = 1.324 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan siswa masuk ke performa yang lebih tinggi, dengan asumsi variabel lain tetap. Dalam parameterisasi polr, ini berasal dari model logit(P(Y <= j)) = alpha_j - x beta.
  • jumlah absensi/ketidakhadiran sekolah memiliki OR kategori lebih tinggi = 0.942 dan p-value = 0.005, sehingga variabel ini menurunkan kecenderungan siswa masuk ke performa yang lebih tinggi, dengan asumsi variabel lain tetap. Dalam parameterisasi polr, ini berasal dari model logit(P(Y <= j)) = alpha_j - x beta.
  • waktu belajar mingguan (skala 1-4) memiliki OR kategori lebih tinggi = 1.291 dan p-value = 0.025, sehingga variabel ini meningkatkan kecenderungan siswa masuk ke performa yang lebih tinggi, dengan asumsi variabel lain tetap. Dalam parameterisasi polr, ini berasal dari model logit(P(Y <= j)) = alpha_j - x beta.
  • jenis kelamin siswa: M dibanding kategori acuan memiliki OR kategori lebih tinggi = 0.640 dan p-value = 0.029, sehingga variabel ini menurunkan kecenderungan siswa masuk ke performa yang lebih tinggi, dengan asumsi variabel lain tetap. Dalam parameterisasi polr, ini berasal dari model logit(P(Y <= j)) = alpha_j - x beta.
cat("Akurasi klasifikasi model ordinal adalah **", scales::percent(akurasi_ordinal),
    "**. AIC model adalah **", fmt_num(AIC_ordinal, 2), "**.\n\n", sep = "")

Akurasi klasifikasi model ordinal adalah 68%. AIC model adalah 1,009.91.

cat("**Kesimpulan substantif:** model ordinal menunjukkan faktor-faktor yang membedakan siswa berkinerja rendah dan tinggi. Variabel seperti jumlah kegagalan kelas sebelumnya, waktu belajar, dan konsumsi alkohol menjadi faktor penting dalam menentukan tingkat capaian akademik. Karena model menggunakan `MASS::polr`, OR di atas 1 berarti variabel tersebut mendorong siswa ke performa yang lebih tinggi.\n")

Kesimpulan substantif: model ordinal menunjukkan faktor-faktor yang membedakan siswa berkinerja rendah dan tinggi. Variabel seperti jumlah kegagalan kelas sebelumnya, waktu belajar, dan konsumsi alkohol menjadi faktor penting dalam menentukan tingkat capaian akademik. Karena model menggunakan MASS::polr, OR di atas 1 berarti variabel tersebut mendorong siswa ke performa yang lebih tinggi.

8 Bagian III — Regresi Logistik Multinomial

8.1 Tujuan Analisis

Tujuan analisis ini adalah mengetahui faktor-faktor yang berhubungan dengan peran jabatan karyawan (JobRole) di perusahaan. Peran jabatan tidak diperlakukan sebagai urutan, karena kategori seperti Sales Executive, Research Scientist, Manager, dan lainnya adalah jenis pekerjaan yang berbeda.

8.2 Variabel yang Digunakan

Variabel respon adalah JobRole:

Variabel yang Dipakai pada Model Multinomial
Kategori JobRole Arti sederhana
Healthcare Representative Perwakilan layanan kesehatan
Human Resources Staf SDM/HR
Laboratory Technician Teknisi laboratorium
Manager Manajer
Manufacturing Director Direktur manufaktur
Research Director Direktur riset
Research Scientist Ilmuwan riset
Sales Executive Eksekutif penjualan (acuan)
Sales Representative Perwakilan penjualan
Variabel Peran Alasan_dipakai
JobRole Respon Mewakili peran jabatan karyawan (multinomial tanpa urutan)
Age Prediktor Usia karyawan berpengaruh pada jenis peran yang diemban
MonthlyIncome Prediktor Pendapatan bulanan mencerminkan level jabatan
TotalWorkingYears Prediktor Total pengalaman kerja mempengaruhi penempatan jabatan
YearsAtCompany Prediktor Lama bekerja di perusahaan
YearsInCurrentRole Prediktor Lama memegang peran jabatan saat ini
Department Prediktor Departemen kerja berkaitan erat dengan peran jabatan
EducationField Prediktor Bidang pendidikan mempengaruhi penempatan jabatan
JobLevel Prediktor Level jabatan mencerminkan hierarki dan peran
Attrition Prediktor Status keluar karyawan mungkin berkaitan dengan peran
BusinessTravel Prediktor Frekuensi perjalanan bisnis bervariasi antar peran
OverTime Prediktor Bekerja lembur dapat bervariasi antar jenis jabatan
Gender Prediktor Jenis kelamin sebagai faktor demografi
MaritalStatus Prediktor Status pernikahan sebagai faktor demografi

Variabel EmployeeNumber, EmployeeCount, Over18, dan StandardHours tidak digunakan karena bersifat konstan atau hanya merupakan identifier. Variabel JobSatisfaction, EnvironmentSatisfaction, WorkLifeBalance, dan JobInvolvement adalah ukuran kepuasan yang mungkin dipengaruhi oleh peran jabatan, bukan sebaliknya, sehingga perlu berhati-hati dalam interpretasi kausalitas.

8.3 Cleaning Data Multinomial

data_multinom <- data_multinom_raw %>%
  dplyr::mutate(
    JobRole        = relevel(factor(JobRole), ref = "Sales Executive"),
    Attrition      = factor(Attrition),
    BusinessTravel = factor(BusinessTravel),
    Department     = factor(Department),
    EducationField = factor(EducationField),
    Gender         = factor(Gender),
    MaritalStatus  = factor(MaritalStatus),
    OverTime       = factor(OverTime)
  ) %>%
  dplyr::select(
    JobRole, Age, MonthlyIncome, TotalWorkingYears, YearsAtCompany,
    YearsInCurrentRole, Department, EducationField, JobLevel,
    Attrition, BusinessTravel, OverTime, Gender, MaritalStatus,
    JobSatisfaction, WorkLifeBalance, EnvironmentSatisfaction,
    DistanceFromHome, PercentSalaryHike, NumCompaniesWorked,
    TrainingTimesLastYear, YearsSinceLastPromotion, YearsWithCurrManager
  ) %>%
  tidyr::drop_na()

knitr::kable(
  head(data_multinom %>%
         dplyr::select(JobRole, Age, MonthlyIncome, Department,
                       EducationField, JobLevel, Attrition, OverTime), 6),
  caption = "Cuplikan Data Multinomial Setelah Cleaning"
)
Cuplikan Data Multinomial Setelah Cleaning
JobRole Age MonthlyIncome Department EducationField JobLevel Attrition OverTime
Sales Executive 41 5993 Sales Life Sciences 2 Yes Yes
Research Scientist 49 5130 Research & Development Life Sciences 2 No No
Laboratory Technician 37 2090 Research & Development Other 1 Yes Yes
Research Scientist 33 2909 Research & Development Life Sciences 1 No Yes
Laboratory Technician 27 3468 Research & Development Medical 1 No No
Laboratory Technician 32 3068 Research & Development Life Sciences 1 No No

Sales Executive digunakan sebagai kategori acuan. Jadi hasil model akan membaca peran jabatan lain dibandingkan dengan Sales Executive.

8.4 Analisis Deskriptif Data Multinomial

dist_multinom <- data_multinom %>%
  dplyr::count(JobRole) %>%
  dplyr::mutate(Persen = scales::percent(n / sum(n))) %>%
  dplyr::arrange(desc(n))
knitr::kable(dist_multinom, caption = "Distribusi Peran Jabatan Karyawan (JobRole)")
Distribusi Peran Jabatan Karyawan (JobRole)
JobRole n Persen
Sales Executive 326 22.18%
Research Scientist 292 19.86%
Laboratory Technician 259 17.62%
Manufacturing Director 145 9.86%
Healthcare Representative 131 8.91%
Manager 102 6.94%
Sales Representative 83 5.65%
Research Director 80 5.44%
Human Resources 52 3.54%
ringkasan_multinom <- data_multinom %>%
  dplyr::summarise(
    Rata_rata_Age          = mean(Age),
    Median_Age             = median(Age),
    Rata_rata_Income       = mean(MonthlyIncome),
    Median_Income          = median(MonthlyIncome),
    Rata_rata_WorkingYears = mean(TotalWorkingYears)
  )
knitr::kable(ringkasan_multinom, caption = "Ringkasan Prediktor Utama Data Multinomial")
Ringkasan Prediktor Utama Data Multinomial
Rata_rata_Age Median_Age Rata_rata_Income Median_Income Rata_rata_WorkingYears
36.92381 36 6502.931 4919 11.27959
ggplot(data_multinom, aes(x = forcats::fct_infreq(JobRole))) +
  geom_bar(fill = "#2f7f73") +
  coord_flip() +
  labs(
    title = "Distribusi Peran Jabatan Karyawan (JobRole)",
    x     = "Job Role",
    y     = "Jumlah karyawan"
  ) +
  theme_minimal()

8.5 Pembentukan Model Regresi Logistik Multinomial

Model multinomial membandingkan setiap kategori JobRole terhadap kategori acuan Sales Executive. Model membangun persamaan logit untuk setiap kategori sekaligus.

Untuk setiap kategori jabatan \(k\) selain Sales Executive, bentuk modelnya adalah:

\[ \log\left( \frac{P(Y_i = k)}{P(Y_i = \text{Sales Executive})} \right) = \beta_{0k} + \beta_{1k}X_{1i} + \cdots + \beta_{pk}X_{pi} \]

model_multinom <- nnet::multinom(
  JobRole ~ Age + MonthlyIncome + TotalWorkingYears + YearsAtCompany +
    YearsInCurrentRole + Department + EducationField + JobLevel +
    Attrition + BusinessTravel + OverTime + Gender + MaritalStatus +
    JobSatisfaction + WorkLifeBalance + EnvironmentSatisfaction +
    DistanceFromHome + PercentSalaryHike + NumCompaniesWorked +
    TrainingTimesLastYear + YearsSinceLastPromotion + YearsWithCurrManager,
  data   = data_multinom,
  trace  = FALSE,
  maxit  = 300,
  MaxNWts= 15000
)

summary(model_multinom)
## Call:
## nnet::multinom(formula = JobRole ~ Age + MonthlyIncome + TotalWorkingYears + 
##     YearsAtCompany + YearsInCurrentRole + Department + EducationField + 
##     JobLevel + Attrition + BusinessTravel + OverTime + Gender + 
##     MaritalStatus + JobSatisfaction + WorkLifeBalance + EnvironmentSatisfaction + 
##     DistanceFromHome + PercentSalaryHike + NumCompaniesWorked + 
##     TrainingTimesLastYear + YearsSinceLastPromotion + YearsWithCurrManager, 
##     data = data_multinom, trace = FALSE, maxit = 300, MaxNWts = 15000)
## 
## Coefficients:
##                           (Intercept)         Age MonthlyIncome
## Healthcare Representative   -193.7347   8.9739256    -0.1375353
## Human Resources             1234.7691   2.6355715    -0.1999099
## Laboratory Technician       -485.0004   8.9783870    -0.1385331
## Manager                     -656.0151 -13.5323613     0.3189753
## Manufacturing Director      -276.9259   8.9579585    -0.1376054
## Research Director          -1170.9908 -13.5773081     0.3188088
## Research Scientist          -286.0277   8.9775036    -0.1385099
## Sales Representative        1032.0450   0.1515729    -0.2168137
##                           TotalWorkingYears YearsAtCompany YearsInCurrentRole
## Healthcare Representative          2.459191       7.424926         -10.194318
## Human Resources                    3.452041     -11.170283          17.246048
## Laboratory Technician              2.577644       7.519987         -10.094700
## Manager                          -11.984791     -12.395691          -3.372370
## Manufacturing Director             2.451742       7.359127         -10.130969
## Research Director                -11.928672     -12.435679          -3.257568
## Research Scientist                 2.586693       7.544561         -10.096615
## Sales Representative              36.454061      -3.872958         -23.728702
##                           DepartmentResearch & Development DepartmentSales
## Healthcare Representative                        795.29932       -430.0472
## Human Resources                                 -891.43044      -2018.5754
## Laboratory Technician                           1137.37243      -1001.6746
## Manager                                         -743.23921      -1165.0631
## Manufacturing Director                           915.03478       -463.8049
## Research Director                                -68.55588       -553.8204
## Research Scientist                              1001.24888      -1434.9438
## Sales Representative                             186.13608        336.0655
##                           EducationFieldLife Sciences EducationFieldMarketing
## Healthcare Representative                  -70.699573              -123.41860
## Human Resources                              5.413813               209.81586
## Laboratory Technician                     -112.121333              -261.14018
## Manager                                   -152.376589                34.05217
## Manufacturing Director                    -107.988059                28.97191
## Research Director                         -308.389921              -224.70220
## Research Scientist                        -174.381050              -271.15342
## Sales Representative                       -28.628892              -127.10677
##                           EducationFieldMedical EducationFieldOther
## Healthcare Representative             -68.16040           -50.99249
## Human Resources                       -61.60781          -204.35700
## Laboratory Technician                -109.97361           -91.98464
## Manager                              -111.96348           -95.97247
## Manufacturing Director               -105.44327           -88.85988
## Research Director                    -268.21754          -253.59824
## Research Scientist                   -172.29094          -154.55773
## Sales Representative                   61.02722          -112.47609
##                           EducationFieldTechnical Degree  JobLevel AttritionYes
## Healthcare Representative                      264.86632  250.3995    -58.47558
## Human Resources                                251.35120  351.2665   -163.20997
## Laboratory Technician                          222.86984  246.2751    -57.05670
## Manager                                         81.70314 -354.3686   -269.12643
## Manufacturing Director                         227.37777  250.8466    -58.19880
## Research Director                              -73.69057 -354.3405   -271.02778
## Research Scientist                             161.48105  245.8634    -58.11508
## Sales Representative                           590.51387 -573.2649   -257.84277
##                           BusinessTravelTravel_Frequently
## Healthcare Representative                        31.66026
## Human Resources                                  51.25161
## Laboratory Technician                            32.34446
## Manager                                        -193.06969
## Manufacturing Director                           32.06820
## Research Director                              -193.13917
## Research Scientist                               32.45998
## Sales Representative                            193.24043
##                           BusinessTravelTravel_Rarely OverTimeYes GenderMale
## Healthcare Representative                    78.97647   -22.28982  -31.89363
## Human Resources                              66.81188    11.24389  -49.84380
## Laboratory Technician                        79.25857   -22.45033  -31.58060
## Manager                                     -99.05467   -66.01441 -107.02564
## Manufacturing Director                       79.30070   -22.31203  -32.36521
## Research Director                           -98.08181   -65.09154 -106.29600
## Research Scientist                           79.48062   -21.66776  -31.87298
## Sales Representative                        329.33854  -116.30310  -16.47816
##                           MaritalStatusMarried MaritalStatusSingle
## Healthcare Representative            224.92532           191.99427
## Human Resources                      -10.40213            80.91409
## Laboratory Technician                225.16249           191.98459
## Manager                               74.28209            39.81838
## Manufacturing Director               224.99214           191.99838
## Research Director                     73.90523            39.80470
## Research Scientist                   225.03891           192.22898
## Sales Representative                  75.43468           184.90523
##                           JobSatisfaction WorkLifeBalance
## Healthcare Representative       -13.34660        59.31338
## Human Resources                 -49.94713       137.59045
## Laboratory Technician           -13.39397        59.59713
## Manager                          87.82190       -14.02239
## Manufacturing Director          -13.39672        59.47529
## Research Director                87.83321       -14.02890
## Research Scientist              -13.32346        59.50167
## Sales Representative             33.94052       -32.11259
##                           EnvironmentSatisfaction DistanceFromHome
## Healthcare Representative               -54.55074       -7.0213057
## Human Resources                         -58.08572       -4.3066132
## Laboratory Technician                   -54.57459       -7.0190072
## Manager                                 -83.74942       -3.1874636
## Manufacturing Director                  -54.39896       -7.0271722
## Research Director                       -84.06256       -3.1434318
## Research Scientist                      -54.64904       -7.0238500
## Sales Representative                    -65.51621        0.3767013
##                           PercentSalaryHike NumCompaniesWorked
## Healthcare Representative          18.96011         -11.061296
## Human Resources                    10.82848         -29.396468
## Laboratory Technician              18.93470         -10.932305
## Manager                            15.68289          -6.823188
## Manufacturing Director             18.97728         -11.072696
## Research Director                  15.66871          -6.623813
## Research Scientist                 18.97797         -10.949900
## Sales Representative               17.49796         -24.087189
##                           TrainingTimesLastYear YearsSinceLastPromotion
## Healthcare Representative            -19.736312               -7.350464
## Human Resources                      -14.094838                7.578287
## Laboratory Technician                -19.646913               -7.453880
## Manager                                5.955875               22.271297
## Manufacturing Director               -19.732722               -7.447505
## Research Director                      5.695241               22.118682
## Research Scientist                   -19.827761               -7.416186
## Sales Representative                  -6.406207               -0.646007
##                           YearsWithCurrManager
## Healthcare Representative            -28.76618
## Human Resources                      -13.22622
## Laboratory Technician                -28.76197
## Manager                              -24.60483
## Manufacturing Director               -28.66626
## Research Director                    -24.54373
## Research Scientist                   -28.80312
## Sales Representative                 -19.70450
## 
## Std. Errors:
##                                 (Intercept)             Age MonthlyIncome
## Healthcare Representative 0.001636025530333 0.0119006104465 0.00007059469
## Human Resources                         NaN             NaN           NaN
## Laboratory Technician     0.002070922882813 0.0105008697927 0.00007596846
## Manager                   0.000587173384428 0.0183056492320 0.00004132232
## Manufacturing Director    0.001700838859934 0.0117261630003 0.00006914167
## Research Director         0.000587173311970 0.0183056440108 0.00004147625
## Research Scientist        0.002102577166302 0.0105565178262 0.00007585756
## Sales Representative      0.000000004943403 0.0000002323387 0.00001471157
##                                   TotalWorkingYears      YearsAtCompany
## Healthcare Representative 0.02387362661423726023502 0.02537669001318218
## Human Resources           0.00000000000000002074766                 NaN
## Laboratory Technician     0.02418860892282232583073 0.03122426517284106
## Manager                   0.02472947126574341278338 0.02276377874454226
## Manufacturing Director    0.02382080419751222283287 0.02493787863774681
## Research Director         0.02472947256476230870614 0.02276377839226153
## Research Scientist        0.02393260740736621178248 0.03074529372671153
## Sales Representative      0.00000002471866503520047 0.00000000001675842
##                                  YearsInCurrentRole
## Healthcare Representative 0.03052644707925256453707
## Human Resources           0.00000000000000001865721
## Laboratory Technician     0.03156598688980761374667
## Manager                   0.03357632803474904642504
## Manufacturing Director    0.03015770760929062030420
## Research Director         0.03357632818128254009737
## Research Scientist        0.03126323567955484328262
## Sales Representative      0.00000000000651716210437
##                           DepartmentResearch & Development
## Healthcare Representative             0.001636029183405088
## Human Resources                                        NaN
## Laboratory Technician                 0.002070922882813129
## Manager                               0.000587173311970338
## Manufacturing Director                0.001700838859483119
## Research Director                     0.000587173311970375
## Research Scientist                    0.002102577166302162
## Sales Representative                  0.000000000000931023
##                                                 DepartmentSales
## Healthcare Representative 0.00000000638827162086904348390833941
## Human Resources                                             NaN
## Laboratory Technician     0.00000000000000000000020519128635730
## Manager                   0.00000000010732756995771270353088461
## Manufacturing Director    0.00000000000263920206070341407727717
## Research Director         0.00000000000000000000000000001813052
## Research Scientist        0.00000000000000000000000000000000000
## Sales Representative      0.00000000494327645416542970930456891
##                           EducationFieldLife Sciences
## Healthcare Representative       0.0008845428188583855
## Human Resources                                   NaN
## Laboratory Technician           0.0025405656761900621
## Manager                         0.0012534732492635056
## Manufacturing Director          0.0010222772619373137
## Research Director               0.0012534732492634698
## Research Scientist              0.0026764714699485884
## Sales Representative            0.0000000000009310232
##                                                                                      EducationFieldMarketing
## Healthcare Representative 0.00000000001094151366904666682396823906842087126278784126043319702148437500000000
## Human Resources           0.00000000000000000477413122486281363537202748759114001586567610502243041992187500
## Laboratory Technician                                                                                    NaN
## Manager                   0.00000000000000000000001241511571628525825109617963803998463845346122980117797852
## Manufacturing Director    0.00000000000263920206070341407727716798170547463087132200598716735839843750000000
## Research Director         0.00000000000000000000000000000000000000000000000000000000000000000000000003822273
## Research Scientist        0.00000000000000000000000000000000000000000000000000000000000000000000000000000000
## Sales Representative      0.00000000000000000000265061710532166948601009992536603476764867082238197326660156
##                           EducationFieldMedical
## Healthcare Representative     0.001560525154756
## Human Resources                             NaN
## Laboratory Technician         0.001867098412612
## Manager                       0.001817337771521
## Manufacturing Director        0.001658248709124
## Research Director             0.001817337771521
## Research Scientist            0.001904897865753
## Sales Representative          0.000000004943276
##                                                                                                       EducationFieldOther
## Healthcare Representative 0.000504691717260775983902587604745804128469899296760559082031250000000000000000000000000000000
## Human Resources                                                                                                       NaN
## Laboratory Technician     0.001063163402480520047799217309147934429347515106201171875000000000000000000000000000000000000
## Manager                   0.000742123304667742229694116229410383311915211379528045654296875000000000000000000000000000000
## Manufacturing Director    0.000415907051641658371231302648851624326198361814022064208984375000000000000000000000000000000
## Research Director         0.000742123304667742446534550726511270113405771553516387939453125000000000000000000000000000000
## Research Scientist        0.001070603243568866297374952978316287044435739517211914062500000000000000000000000000000000000
## Sales Representative      0.000000000000000000000000000000000000000000000000000000000000000000000000000000000000002517525
##                           EducationFieldTechnical Degree
## Healthcare Representative   0.00064060429485782621923662
## Human Resources                                      NaN
## Laboratory Technician       0.00134278618064070085803829
## Manager                     0.00092682365934712357055658
## Manufacturing Director      0.00059318261718578808600749
## Research Director           0.00092682360980200762228670
## Research Scientist          0.00144881331639602740045392
## Sales Representative        0.00000000000000000004873736
##                                            JobLevel
## Healthcare Representative 0.00205699073121571500769
## Human Resources           0.00000000000000002914509
## Laboratory Technician     0.00151841061368537789088
## Manager                   0.00070455836363053875142
## Manufacturing Director    0.00215055573440868248031
## Research Director         0.00070455824250549830678
## Research Scientist        0.00162570040323980363230
## Sales Representative      0.00000000988667945900040
##                                           AttritionYes
## Healthcare Representative 0.00075810334407961472383969
## Human Resources                                    NaN
## Laboratory Technician     0.00199288243452218610113347
## Manager                   0.00044450665026605817907804
## Manufacturing Director    0.00075284637701391839744586
## Research Director         0.00044450665026605817907804
## Research Scientist        0.00176686046629414617893661
## Sales Representative      0.00000000000000000000276749
##                           BusinessTravelTravel_Frequently
## Healthcare Representative               0.001016306935012
## Human Resources                                       NaN
## Laboratory Technician                   0.001493608623050
## Manager                                 0.000916612847806
## Manufacturing Director                  0.000874829397026
## Research Director                       0.000916612847806
## Research Scientist                      0.001410495145240
## Sales Representative                    0.000000004943276
##                           BusinessTravelTravel_Rarely
## Healthcare Representative  0.001628706137962928788762
## Human Resources            0.000000000000000001088456
## Laboratory Technician      0.002052301019080409860001
## Manager                    0.000568828934115224145798
## Manufacturing Director     0.001635578296240736313749
## Research Director          0.000568828850159129898942
## Research Scientist         0.002004366972167395243715
## Sales Representative       0.000000000000931023052908
##                                              OverTimeYes
## Healthcare Representative 0.0006797489045432548440336951
## Human Resources                                      NaN
## Laboratory Technician     0.0021428841740646622712751501
## Manager                   0.0019207375782153575903554721
## Manufacturing Director    0.0005678841930602837421993012
## Research Director         0.0019207374879906787337735530
## Research Scientist        0.0022334226167275263510991756
## Sales Representative      0.0000000000000000000003393096
##                                         GenderMale
## Healthcare Representative 0.0010257715808876994613
## Human Resources           0.0000000000000000477711
## Laboratory Technician     0.0041382449293419498737
## Manager                   0.0012302687527956453380
## Manufacturing Director    0.0008377975329476091554
## Research Director         0.0012302687573596303786
## Research Scientist        0.0041032135366639117466
## Sales Representative      0.0000000049432764535056
##                                     MaritalStatusMarried
## Healthcare Representative 0.0015205226581476033193701669
## Human Resources                                      NaN
## Laboratory Technician     0.0016468321399982083219626361
## Manager                   0.0018976766897374022259603610
## Manufacturing Director    0.0015410703304662427652604872
## Research Director         0.0018976766927394060364281403
## Research Scientist        0.0019902173909439784330044354
## Sales Representative      0.0000000000000000000003393096
##                                  MaritalStatusSingle  JobSatisfaction
## Healthcare Representative 0.001411449394486838761556 0.00893059191295
## Human Resources           0.000000000000000009002845              NaN
## Laboratory Technician     0.002684128858123998109808 0.03791829075527
## Manager                   0.001200400435271526625863 0.00504277495354
## Manufacturing Director    0.001664901327257130780532 0.00874554171666
## Research Director         0.001200400435271549177269 0.00504277486374
## Research Scientist        0.002861877633423594909784 0.03782582266990
## Sales Representative      0.000000004943276453077128 0.00000001482996
##                            WorkLifeBalance EnvironmentSatisfaction
## Healthcare Representative 0.00547548034186        0.01139910817194
## Human Resources                        NaN                     NaN
## Laboratory Technician     0.00845676629116        0.03699610529655
## Manager                   0.00224101470282        0.00207440467826
## Manufacturing Director    0.00571756094418        0.01166427156443
## Research Director         0.00224101470289        0.00207440478203
## Research Scientist        0.00862005858322        0.03700484664013
## Sales Representative      0.00000001483008        0.00000000988706
##                                      DistanceFromHome
## Healthcare Representative 0.0110799574626983237529698
## Human Resources           0.0000000000000000009644243
## Laboratory Technician     0.0098012203090086895407307
## Manager                   0.0139714346613242224037776
## Manufacturing Director    0.0109086601371174169577971
## Research Director         0.0139714346914975948849014
## Research Scientist        0.0098925655351133935655472
## Sales Representative      0.0000000197748811811796246
##                                    PercentSalaryHike NumCompaniesWorked
## Healthcare Representative 0.022842101846402788206669   0.04092498436839
## Human Resources           0.000000000000000007724598                NaN
## Laboratory Technician     0.019848978120084289494818   0.03579188579651
## Manager                   0.023626888478305234364552   0.04508776702955
## Manufacturing Director    0.021688991152221732960914   0.03993711255418
## Research Director         0.023626888821036241133955   0.04508776708183
## Research Scientist        0.019642981068175021514000   0.03582250349854
## Sales Representative      0.000000093923771402119989   0.00000001482996
##                           TrainingTimesLastYear YearsSinceLastPromotion
## Healthcare Representative      0.03941278442076     0.03856968680186983
## Human Resources                             NaN                     NaN
## Laboratory Technician          0.03289107784399     0.03863004946199855
## Manager                        0.00460650650370     0.03232816839753523
## Manufacturing Director         0.04056147162006     0.03898866456054612
## Research Director              0.00460650644860     0.03232816839753531
## Research Scientist             0.03238756166211     0.03820453001122903
## Sales Representative           0.00000001483008     0.00000000001117228
##                           YearsWithCurrManager
## Healthcare Representative  0.03144712159051799
## Human Resources                            NaN
## Laboratory Technician      0.03141118633364032
## Manager                    0.03645186601656455
## Manufacturing Director     0.03074123837897624
## Research Director          0.03645186611473998
## Research Scientist         0.03172568501362495
## Sales Representative       0.00000000001582739
## 
## Residual Deviance: 1541.73 
## AIC: 2021.73

8.6 Persamaan Model Hasil Estimasi Multinomial

Chunk model-multinomial di atas sudah membangun model regresi logistik multinomial menggunakan nnet::multinom(). Persamaan logit yang dibentuk adalah untuk setiap kategori JobRole dibandingkan Sales Executive, misalnya:

\[ \log\left(\frac{P(Y=\mathrm{Manager})}{P(Y=\mathrm{Sales\ Executive})}\right) = \beta_{0,\mathrm{Manager}} + \beta_{1,\mathrm{Manager}}\,\mathrm{Age} + \cdots + \beta_{p,\mathrm{Manager}}\,X_p \]

8.7 Koefisien dan Relative Risk Ratio Model Multinomial

summary_multinom <- summary(model_multinom)
coef_multinom    <- summary_multinom$coefficients
se_multinom      <- summary_multinom$standard.errors

multinom_table <- purrr::map_dfr(rownames(coef_multinom), function(kelas) {
  tibble::tibble(
    Kategori  = kelas,
    term      = colnames(coef_multinom),
    estimate  = as.numeric(coef_multinom[kelas, ]),
    std_error = as.numeric(se_multinom[kelas, ])
  )
}) %>%
  dplyr::mutate(
    z_value              = estimate / std_error,
    p_value              = 2 * pnorm(abs(z_value), lower.tail = FALSE),
    Relative_Risk_Ratio  = exp(estimate),
    Variabel_Mudah       = purrr::map_chr(term, label_multinom),
    Keputusan_5persen    = ifelse(p_value < 0.05, "Signifikan", "Tidak signifikan")
  ) %>%
  dplyr::arrange(Kategori, p_value)

knitr::kable(
  multinom_table %>%
    dplyr::select(Kategori, Variabel_Mudah, term, estimate, std_error,
                  z_value, p_value, Relative_Risk_Ratio, Keputusan_5persen) %>%
    head(40),
  digits  = 4,
  caption = "Koefisien dan Relative Risk Ratio Regresi Multinomial"
)
Koefisien dan Relative Risk Ratio Regresi Multinomial
Kategori Variabel_Mudah term estimate std_error z_value p_value Relative_Risk_Ratio Keputusan_5persen
Healthcare Representative (Intercept) (Intercept) -193.7347 0.0016 -118417.8912 0 0.0000 Signifikan
Healthcare Representative usia karyawan Age 8.9739 0.0119 754.0727 0 7894.5316 Signifikan
Healthcare Representative pendapatan bulanan karyawan MonthlyIncome -0.1375 0.0001 -1948.2385 0 0.8715 Signifikan
Healthcare Representative total masa kerja (tahun) TotalWorkingYears 2.4592 0.0239 103.0087 0 11.6954 Signifikan
Healthcare Representative lama bekerja di perusahaan saat ini YearsAtCompany 7.4249 0.0254 292.5884 0 1677.2749 Signifikan
Healthcare Representative lama memegang peran jabatan saat ini YearsInCurrentRole -10.1943 0.0305 -333.9504 0 0.0000 Signifikan
Healthcare Representative departemen kerja: Research & Development dibanding acuan DepartmentResearch & Development 795.2993 0.0016 486115.6060 0 Inf Signifikan
Healthcare Representative departemen kerja: Sales dibanding acuan DepartmentSales -430.0472 0.0000 -67318234280.1656 0 0.0000 Signifikan
Healthcare Representative bidang pendidikan: Life Sciences dibanding acuan EducationFieldLife Sciences -70.6996 0.0009 -79927.8130 0 0.0000 Signifikan
Healthcare Representative bidang pendidikan: Marketing dibanding acuan EducationFieldMarketing -123.4186 0.0000 -11279847252144.3848 0 0.0000 Signifikan
Healthcare Representative bidang pendidikan: Medical dibanding acuan EducationFieldMedical -68.1604 0.0016 -43677.8606 0 0.0000 Signifikan
Healthcare Representative bidang pendidikan: Other dibanding acuan EducationFieldOther -50.9925 0.0005 -101036.9016 0 0.0000 Signifikan
Healthcare Representative bidang pendidikan: Technical Degree dibanding acuan EducationFieldTechnical Degree 264.8663 0.0006 413463.2277 0 10714723834898563653422080068024864426208660466464480806626082860882400662820204660080062440040420860606468484480606.0000 Signifikan
Healthcare Representative level jabatan karyawan JobLevel 250.3995 0.0021 121730.9990 0 5586464935741923457206606602042888460440244608880482264006488082844064864246288682468880042606224002004862202.0000 Signifikan
Healthcare Representative status keluar karyawan: Yes dibanding acuan AttritionYes -58.4756 0.0008 -77134.0468 0 0.0000 Signifikan
Healthcare Representative frekuensi perjalanan bisnis: Travel_Frequently dibanding acuan BusinessTravelTravel_Frequently 31.6603 0.0010 31152.2620 0 56218101026456.5156 Signifikan
Healthcare Representative frekuensi perjalanan bisnis: Travel_Rarely dibanding acuan BusinessTravelTravel_Rarely 78.9765 0.0016 48490.3152 0 19908882855787660328482624008408666.0000 Signifikan
Healthcare Representative bekerja lembur: Yes dibanding acuan OverTimeYes -22.2898 0.0007 -32791.2601 0 0.0000 Signifikan
Healthcare Representative jenis kelamin karyawan: Male dibanding acuan GenderMale -31.8936 0.0010 -31092.3280 0 0.0000 Signifikan
Healthcare Representative status pernikahan: Married dibanding acuan MaritalStatusMarried 224.9253 0.0015 147926.3188 0 48286430097978824912026266200400440486828886444864062862208846680668886802408260848040248406084860.0000 Signifikan
Healthcare Representative status pernikahan: Single dibanding acuan MaritalStatusSingle 191.9943 0.0014 136026.3191 0 241018371863804662582200004664686062268446842024062804624648808426848004600880282084.0000 Signifikan
Healthcare Representative tingkat kepuasan kerja (skala 1-4) JobSatisfaction -13.3466 0.0089 -1494.4811 0 0.0000 Signifikan
Healthcare Representative keseimbangan hidup-kerja (skala 1-4) WorkLifeBalance 59.3134 0.0055 10832.5441 0 57474425417065050072440480.0000 Signifikan
Healthcare Representative kepuasan lingkungan kerja (skala 1-4) EnvironmentSatisfaction -54.5507 0.0114 -4785.5269 0 0.0000 Signifikan
Healthcare Representative jarak dari rumah ke kantor (km) DistanceFromHome -7.0213 0.0111 -633.6943 0 0.0009 Signifikan
Healthcare Representative persentase kenaikan gaji PercentSalaryHike 18.9601 0.0228 830.0509 0 171502285.1237 Signifikan
Healthcare Representative jumlah perusahaan sebelumnya NumCompaniesWorked -11.0613 0.0409 -270.2822 0 0.0000 Signifikan
Healthcare Representative frekuensi pelatihan tahun lalu TrainingTimesLastYear -19.7363 0.0394 -500.7592 0 0.0000 Signifikan
Healthcare Representative tahun sejak promosi terakhir YearsSinceLastPromotion -7.3505 0.0386 -190.5762 0 0.0006 Signifikan
Healthcare Representative lama dengan manajer saat ini YearsWithCurrManager -28.7662 0.0314 -914.7478 0 0.0000 Signifikan
Human Resources total masa kerja (tahun) TotalWorkingYears 3.4520 0.0000 166382212338542496.0000 0 31.5648 Signifikan
Human Resources lama memegang peran jabatan saat ini YearsInCurrentRole 17.2460 0.0000 924363539835938176.0000 0 30893240.0541 Signifikan
Human Resources bidang pendidikan: Marketing dibanding acuan EducationFieldMarketing 209.8159 0.0000 43948490017981497344.0000 0 13239438112204088689824262462024004266822026442288022684220462664422686806864624240464426282.0000 Signifikan
Human Resources level jabatan karyawan JobLevel 351.2665 0.0000 12052341435678699520.0000 0 357368189835990547028220622822246644268288082062806224040860606264640640042686446800226604084864026648284442686446428000060040028028806282624664840648288.0000 Signifikan
Human Resources frekuensi perjalanan bisnis: Travel_Rarely dibanding acuan BusinessTravelTravel_Rarely 66.8119 0.0000 61382254040884273152.0000 0 103759853300141322560488404226.0000 Signifikan
Human Resources jenis kelamin karyawan: Male dibanding acuan GenderMale -49.8438 0.0000 -1043388109321112064.0000 0 0.0000 Signifikan
Human Resources status pernikahan: Single dibanding acuan MaritalStatusSingle 80.9141 0.0000 8987613563558907904.0000 0 138210813838098451120288882046268028.0000 Signifikan
Human Resources jarak dari rumah ke kantor (km) DistanceFromHome -4.3066 0.0000 -4465475619272303616.0000 0 0.0135 Signifikan
Human Resources persentase kenaikan gaji PercentSalaryHike 10.8285 0.0000 1401818318499469312.0000 0 50437.1094 Signifikan
Human Resources (Intercept) (Intercept) 1234.7691 NaN NaN NaN Inf NA

8.8 Evaluasi Model Multinomial

pred_multinom <- predict(model_multinom, type = "class")
conf_multinom <- table(Prediksi = pred_multinom, Aktual = data_multinom$JobRole)
conf_multinom
##                            Aktual
## Prediksi                    Sales Executive Healthcare Representative
##   Sales Executive                       326                         0
##   Healthcare Representative               0                        65
##   Human Resources                         0                         0
##   Laboratory Technician                   0                         8
##   Manager                                 0                         0
##   Manufacturing Director                  0                        51
##   Research Director                       0                         0
##   Research Scientist                      0                         7
##   Sales Representative                    0                         0
##                            Aktual
## Prediksi                    Human Resources Laboratory Technician Manager
##   Sales Executive                         0                     0       0
##   Healthcare Representative               0                    12       0
##   Human Resources                        52                     0       0
##   Laboratory Technician                   0                   132       0
##   Manager                                 0                     0      78
##   Manufacturing Director                  0                    10       0
##   Research Director                       0                     0      24
##   Research Scientist                      0                   105       0
##   Sales Representative                    0                     0       0
##                            Aktual
## Prediksi                    Manufacturing Director Research Director
##   Sales Executive                                0                 0
##   Healthcare Representative                     37                 0
##   Human Resources                                0                 0
##   Laboratory Technician                          9                 0
##   Manager                                        0                14
##   Manufacturing Director                        89                 0
##   Research Director                              0                66
##   Research Scientist                            10                 0
##   Sales Representative                           0                 0
##                            Aktual
## Prediksi                    Research Scientist Sales Representative
##   Sales Executive                            0                    0
##   Healthcare Representative                  5                    0
##   Human Resources                            0                    0
##   Laboratory Technician                     75                    0
##   Manager                                    0                    0
##   Manufacturing Director                    15                    0
##   Research Director                          0                    0
##   Research Scientist                       197                    0
##   Sales Representative                       0                   83
akurasi_multinom <- mean(as.character(pred_multinom) == as.character(data_multinom$JobRole))
akurasi_multinom
## [1] 0.7401361
AIC_multinom <- AIC(model_multinom)
AIC_multinom
## [1] 2021.73
prob_multinom <- predict(model_multinom, type = "probs")
knitr::kable(
  head(prob_multinom, 6),
  digits  = 4,
  caption = "Cuplikan Prediksi Probabilitas untuk Model Multinomial"
)
Cuplikan Prediksi Probabilitas untuk Model Multinomial
Sales Executive Healthcare Representative Human Resources Laboratory Technician Manager Manufacturing Director Research Director Research Scientist Sales Representative
1 0.0000 0 0.0000 0 0.0000 0 0.0000 0
0 0.1600 0 0.3373 0 0.2671 0 0.2357 0
0 0.0003 0 0.6314 0 0.0002 0 0.3681 0
0 0.0007 0 0.2911 0 0.0007 0 0.7075 0
0 0.0034 0 0.5798 0 0.0022 0 0.4145 0
0 0.0017 0 0.4071 0 0.0019 0 0.5893 0

8.9 Visualisasi Probabilitas Prediksi Model Multinomial

Visualisasi probabilitas prediksi digunakan untuk melihat bagaimana peluang masing-masing peran jabatan berubah ketika MonthlyIncome berubah. Variabel lain ditahan pada nilai yang umum.

grid_multinom <- data.frame(
  MonthlyIncome = seq(
    quantile(data_multinom$MonthlyIncome, 0.05, na.rm = TRUE),
    quantile(data_multinom$MonthlyIncome, 0.95, na.rm = TRUE),
    length.out  = 100
  ),
  Age                    = median(data_multinom$Age, na.rm = TRUE),
  TotalWorkingYears      = median(data_multinom$TotalWorkingYears, na.rm = TRUE),
  YearsAtCompany         = median(data_multinom$YearsAtCompany, na.rm = TRUE),
  YearsInCurrentRole     = median(data_multinom$YearsInCurrentRole, na.rm = TRUE),
  Department             = factor(modus(data_multinom$Department),    levels = levels(data_multinom$Department)),
  EducationField         = factor(modus(data_multinom$EducationField), levels = levels(data_multinom$EducationField)),
  JobLevel               = median(data_multinom$JobLevel, na.rm = TRUE),
  Attrition              = factor(modus(data_multinom$Attrition),      levels = levels(data_multinom$Attrition)),
  BusinessTravel         = factor(modus(data_multinom$BusinessTravel), levels = levels(data_multinom$BusinessTravel)),
  OverTime               = factor(modus(data_multinom$OverTime),       levels = levels(data_multinom$OverTime)),
  Gender                 = factor(modus(data_multinom$Gender),         levels = levels(data_multinom$Gender)),
  MaritalStatus          = factor(modus(data_multinom$MaritalStatus),  levels = levels(data_multinom$MaritalStatus)),
  JobSatisfaction        = median(data_multinom$JobSatisfaction, na.rm = TRUE),
  WorkLifeBalance        = median(data_multinom$WorkLifeBalance, na.rm = TRUE),
  EnvironmentSatisfaction= median(data_multinom$EnvironmentSatisfaction, na.rm = TRUE),
  DistanceFromHome       = median(data_multinom$DistanceFromHome, na.rm = TRUE),
  PercentSalaryHike      = median(data_multinom$PercentSalaryHike, na.rm = TRUE),
  NumCompaniesWorked     = median(data_multinom$NumCompaniesWorked, na.rm = TRUE),
  TrainingTimesLastYear  = median(data_multinom$TrainingTimesLastYear, na.rm = TRUE),
  YearsSinceLastPromotion= median(data_multinom$YearsSinceLastPromotion, na.rm = TRUE),
  YearsWithCurrManager   = median(data_multinom$YearsWithCurrManager, na.rm = TRUE)
)

prob_grid_multinom <- predict(model_multinom, newdata = grid_multinom, type = "probs")

plot_prob_multinom <- grid_multinom %>%
  dplyr::bind_cols(as.data.frame(prob_grid_multinom)) %>%
  tidyr::pivot_longer(
    cols      = levels(data_multinom$JobRole),
    names_to  = "JobRole",
    values_to = "Probabilitas"
  )

ggplot(
  plot_prob_multinom,
  aes(x = MonthlyIncome, y = Probabilitas, color = JobRole)
) +
  geom_line(linewidth = 0.9) +
  scale_y_continuous(labels = scales::percent_format()) +
  labs(
    title    = "Prediksi Probabilitas Peran Jabatan Berdasarkan Pendapatan Bulanan",
    subtitle = "Variabel lain ditahan pada median atau kategori paling umum",
    x        = "Pendapatan bulanan",
    y        = "Probabilitas prediksi",
    color    = "Job Role"
  ) +
  theme_minimal() +
  theme(legend.position = "bottom")

Grafik ini membantu membaca hasil regresi multinomial secara lebih substantif. Jika probabilitas kategori Manager atau Research Director meningkat ketika MonthlyIncome naik, maka model menunjukkan bahwa karyawan dengan pendapatan lebih tinggi cenderung menduduki jabatan manajerial atau direktoral. Sebaliknya, jabatan yang probabilitasnya menurun seiring kenaikan pendapatan cenderung dikaitkan dengan level pendapatan yang lebih rendah.

8.10 Catatan Evaluasi dan Asumsi Model Multinomial

Pada analisis ini, uji goodness-of-fit devians \(G^2\) tidak digunakan sebagai evaluasi utama karena data HR berbentuk data individu dengan beberapa prediktor kontinu. Uji \(G^2\) lebih cocok digunakan pada data berbentuk grouped data atau tabel frekuensi. Evaluasi model multinomial dalam laporan ini difokuskan pada koefisien model, relative risk ratio, AIC, prediksi probabilitas, dan confusion matrix.

Regresi logistik multinomial memiliki asumsi Independence of Irrelevant Alternatives (IIA). Dalam konteks peran jabatan, asumsi ini menyatakan bahwa perbandingan peluang antara dua kategori JobRole tidak berubah secara drastis hanya karena kategori lain ditambahkan atau dihapus. Asumsi IIA dicatat sebagai keterbatasan dalam laporan ini.

8.11 Interpretasi dan Kesimpulan Model Multinomial

sig_multi <- multinom_table %>%
  dplyr::filter(term != "(Intercept)", p_value < 0.05) %>%
  dplyr::arrange(Kategori, p_value)

cat("Pada taraf signifikansi 5%, model multinomial menemukan ", nrow(sig_multi),
    " kombinasi prediktor-kategori yang signifikan. Semua interpretasi dibandingkan terhadap kategori acuan `Sales Executive`.\n\n", sep = "")

Pada taraf signifikansi 5%, model multinomial menemukan 211 kombinasi prediktor-kategori yang signifikan. Semua interpretasi dibandingkan terhadap kategori acuan Sales Executive.

if (nrow(sig_multi) > 0) {
  for (kelas in unique(sig_multi$Kategori)) {
    cat("**JobRole ", kelas, " dibanding Sales Executive**\n\n", sep = "")
    top_kelas <- sig_multi %>% dplyr::filter(Kategori == kelas) %>% head(5)
    for (i in seq_len(nrow(top_kelas))) {
      row  <- top_kelas[i, ]
      arah <- ifelse(row$Relative_Risk_Ratio > 1,
                     paste0("meningkatkan kecenderungan jabatan ", kelas, " dibanding Sales Executive"),
                     paste0("menurunkan kecenderungan jabatan ", kelas, " dibanding Sales Executive"))
      cat("- **", row$Variabel_Mudah, "** memiliki RRR = ", fmt_num(row$Relative_Risk_Ratio, 3),
          " dan p-value = ", fmt_p(row$p_value), ", sehingga variabel ini ", arah,
          ", dengan asumsi variabel lain tetap.\n", sep = "")
    }
    cat("\n")
  }
}

JobRole Healthcare Representative dibanding Sales Executive

  • usia karyawan memiliki RRR = 7,894.532 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Healthcare Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • pendapatan bulanan karyawan memiliki RRR = 0.872 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Healthcare Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • total masa kerja (tahun) memiliki RRR = 11.695 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Healthcare Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama bekerja di perusahaan saat ini memiliki RRR = 1,677.275 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Healthcare Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama memegang peran jabatan saat ini memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Healthcare Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.

JobRole Human Resources dibanding Sales Executive

  • total masa kerja (tahun) memiliki RRR = 31.565 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Human Resources dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama memegang peran jabatan saat ini memiliki RRR = 30,893,240.054 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Human Resources dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • bidang pendidikan: Marketing dibanding acuan memiliki RRR = 13,239,438,112,204,088,689,824,262,462,024,004,266,822,026,442,288,022,684,220,462,664,422,686,806,864,624,240,464,426,282.000 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Human Resources dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • level jabatan karyawan memiliki RRR = 357,368,189,835,990,547,028,220,622,822,246,644,268,288,082,062,806,224,040,860,606,264,640,640,042,686,446,800,226,604,084,864,026,648,284,442,686,446,428,000,060,040,028,028,806,282,624,664,840,648,288.000 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Human Resources dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • frekuensi perjalanan bisnis: Travel_Rarely dibanding acuan memiliki RRR = 103,759,853,300,141,322,560,488,404,226.000 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Human Resources dibanding Sales Executive, dengan asumsi variabel lain tetap.

JobRole Laboratory Technician dibanding Sales Executive

  • usia karyawan memiliki RRR = 7,929.831 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Laboratory Technician dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • pendapatan bulanan karyawan memiliki RRR = 0.871 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Laboratory Technician dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • total masa kerja (tahun) memiliki RRR = 13.166 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Laboratory Technician dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama bekerja di perusahaan saat ini memiliki RRR = 1,844.543 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Laboratory Technician dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama memegang peran jabatan saat ini memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Laboratory Technician dibanding Sales Executive, dengan asumsi variabel lain tetap.

JobRole Manager dibanding Sales Executive

  • usia karyawan memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Manager dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • pendapatan bulanan karyawan memiliki RRR = 1.376 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Manager dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • total masa kerja (tahun) memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Manager dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama bekerja di perusahaan saat ini memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Manager dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama memegang peran jabatan saat ini memiliki RRR = 0.034 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Manager dibanding Sales Executive, dengan asumsi variabel lain tetap.

JobRole Manufacturing Director dibanding Sales Executive

  • usia karyawan memiliki RRR = 7,769.480 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Manufacturing Director dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • pendapatan bulanan karyawan memiliki RRR = 0.871 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Manufacturing Director dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • total masa kerja (tahun) memiliki RRR = 11.609 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Manufacturing Director dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama bekerja di perusahaan saat ini memiliki RRR = 1,570.465 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Manufacturing Director dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama memegang peran jabatan saat ini memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Manufacturing Director dibanding Sales Executive, dengan asumsi variabel lain tetap.

JobRole Research Director dibanding Sales Executive

  • usia karyawan memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Research Director dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • pendapatan bulanan karyawan memiliki RRR = 1.375 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Research Director dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • total masa kerja (tahun) memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Research Director dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama bekerja di perusahaan saat ini memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Research Director dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama memegang peran jabatan saat ini memiliki RRR = 0.038 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Research Director dibanding Sales Executive, dengan asumsi variabel lain tetap.

JobRole Research Scientist dibanding Sales Executive

  • usia karyawan memiliki RRR = 7,922.829 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Research Scientist dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • pendapatan bulanan karyawan memiliki RRR = 0.871 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Research Scientist dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • total masa kerja (tahun) memiliki RRR = 13.286 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Research Scientist dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama bekerja di perusahaan saat ini memiliki RRR = 1,890.432 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Research Scientist dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama memegang peran jabatan saat ini memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Research Scientist dibanding Sales Executive, dengan asumsi variabel lain tetap.

JobRole Sales Representative dibanding Sales Executive

  • usia karyawan memiliki RRR = 1.164 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Sales Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • pendapatan bulanan karyawan memiliki RRR = 0.805 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Sales Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • total masa kerja (tahun) memiliki RRR = 6,788,867,820,348,671.000 dan p-value = < 0.001, sehingga variabel ini meningkatkan kecenderungan jabatan Sales Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama bekerja di perusahaan saat ini memiliki RRR = 0.021 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Sales Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
  • lama memegang peran jabatan saat ini memiliki RRR = 0.000 dan p-value = < 0.001, sehingga variabel ini menurunkan kecenderungan jabatan Sales Representative dibanding Sales Executive, dengan asumsi variabel lain tetap.
cat("Akurasi klasifikasi model multinomial adalah **", scales::percent(akurasi_multinom),
    "**. AIC model adalah **", fmt_num(AIC_multinom, 2), "**.\n\n", sep = "")

Akurasi klasifikasi model multinomial adalah 74%. AIC model adalah 2,021.73.

cat("**Kesimpulan substantif:** model multinomial menunjukkan karakteristik karyawan yang membedakan peran jabatan mereka. Faktor seperti pendapatan bulanan, departemen, total pengalaman kerja, dan bidang pendidikan merupakan pembeda utama antar kategori `JobRole`. Faktor signifikan dengan RRR lebih dari 1 meningkatkan kecenderungan suatu jabatan dibanding `Sales Executive`, sedangkan RRR kurang dari 1 menurunkannya.\n")

Kesimpulan substantif: model multinomial menunjukkan karakteristik karyawan yang membedakan peran jabatan mereka. Faktor seperti pendapatan bulanan, departemen, total pengalaman kerja, dan bidang pendidikan merupakan pembeda utama antar kategori JobRole. Faktor signifikan dengan RRR lebih dari 1 meningkatkan kecenderungan suatu jabatan dibanding Sales Executive, sedangkan RRR kurang dari 1 menurunkannya.

9 Bagian IV — Regresi Poisson

9.1 Tujuan Analisis

Tujuan analisis ini adalah mengetahui faktor-faktor yang berhubungan dengan jumlah dokter yang dikunjungi oleh pasien lansia. Data bersumber dari National Poll on Healthy Aging (NPHA). Variabel respon berupa data cacah (bilangan bulat non-negatif) sehingga Regresi Poisson merupakan pilihan model yang tepat.

9.2 Variabel yang Digunakan

Variabel yang Dipakai pada Model Poisson
Variabel Peran Alasan_dipakai
Number of Doctors Visited Respon Jumlah dokter yang dikunjungi pasien lansia (data cacah)
Age Prediktor Kelompok usia pasien dapat memengaruhi frekuensi kunjungan
Phyiscal Health Prediktor Status kesehatan fisik berkaitan dengan kebutuhan kunjungan
Mental Health Prediktor Status kesehatan mental dapat mendorong atau menghambat kunjungan
Dental Health Prediktor Status kesehatan gigi sebagai dimensi kesehatan lain
Employment Prediktor Status pekerjaan mempengaruhi akses dan waktu kunjungan
Trouble Sleeping Prediktor Gangguan tidur umum mungkin berkaitan dengan jumlah kondisi medis
Stress Keeps Patient from Sleeping Prediktor Gangguan tidur akibat stres
Medication Keeps Patient from Sleeping Prediktor Gangguan tidur akibat efek obat (indikator penggunaan obat)
Pain Keeps Patient from Sleeping Prediktor Gangguan tidur akibat nyeri fisik
Bathroom Needs Keeps Patient from Sleeping Prediktor Gangguan tidur akibat kebutuhan ke toilet
Uknown Keeps Patient from Sleeping Prediktor Gangguan tidur penyebab tidak diketahui
Prescription Sleep Medication Prediktor Penggunaan obat tidur resep (indikator konsumsi obat)
Race Prediktor Kategori ras pasien sebagai faktor sosio-demografis
Gender Prediktor Jenis kelamin pasien

9.3 Cleaning Data Poisson

# Bersihkan nama kolom agar mudah digunakan dalam R
data_poisson <- data_poisson_raw %>%
  dplyr::rename(
    Num_Doctors_Visited                      = `Number of Doctors Visited`,
    Age                                      = Age,
    Physical_Health                          = `Phyiscal Health`,
    Mental_Health                            = `Mental Health`,
    Dental_Health                            = `Dental Health`,
    Employment                               = Employment,
    Stress_Sleep                             = `Stress Keeps Patient from Sleeping`,
    Medication_Sleep                         = `Medication Keeps Patient from Sleeping`,
    Pain_Sleep                               = `Pain Keeps Patient from Sleeping`,
    Bathroom_Sleep                           = `Bathroom Needs Keeps Patient from Sleeping`,
    Unknown_Sleep                            = `Uknown Keeps Patient from Sleeping`,
    Trouble_Sleeping                         = `Trouble Sleeping`,
    Prescription_Sleep_Med                   = `Prescription Sleep Medication`,
    Race                                     = Race,
    Gender                                   = Gender
  ) %>%
  dplyr::mutate(
    Num_Doctors_Visited   = as.integer(Num_Doctors_Visited),
    Age                   = as.integer(Age),
    Physical_Health       = as.integer(Physical_Health),
    Mental_Health         = as.integer(Mental_Health),
    Dental_Health         = as.integer(Dental_Health),
    Employment            = as.integer(Employment),
    Stress_Sleep          = as.integer(Stress_Sleep),
    Medication_Sleep      = as.integer(Medication_Sleep),
    Pain_Sleep            = as.integer(Pain_Sleep),
    Bathroom_Sleep        = as.integer(Bathroom_Sleep),
    Unknown_Sleep         = as.integer(Unknown_Sleep),
    Trouble_Sleeping      = as.integer(Trouble_Sleeping),
    Prescription_Sleep_Med= as.integer(Prescription_Sleep_Med),
    Race                  = as.integer(Race),
    Gender                = as.integer(Gender)
  ) %>%
  tidyr::drop_na()

knitr::kable(
  head(data_poisson, 8),
  caption = "Cuplikan Data Poisson Setelah Cleaning"
)
Cuplikan Data Poisson Setelah Cleaning
Num_Doctors_Visited Age Physical_Health Mental_Health Dental_Health Employment Stress_Sleep Medication_Sleep Pain_Sleep Bathroom_Sleep Unknown_Sleep Trouble_Sleeping Prescription_Sleep_Med Race Gender
3 2 4 3 3 3 0 0 0 0 1 2 3 1 2
2 2 4 2 3 3 1 0 0 1 0 3 3 1 1
3 2 3 2 3 3 0 0 0 0 1 3 3 4 1
1 2 3 2 3 3 0 0 0 1 0 3 3 4 2
3 2 3 3 3 3 1 0 0 0 0 2 3 1 2
2 2 3 2 4 3 0 0 0 1 0 3 3 1 1
3 2 4 1 1 3 0 0 1 1 0 2 1 1 1
2 2 3 2 6 3 1 0 0 0 0 3 3 1 1

9.4 Analisis Deskriptif Data Poisson

ringkasan_poisson <- data_poisson %>%
  dplyr::summarise(
    Jumlah_observasi      = dplyr::n(),
    Rata_rata_kunjungan   = mean(Num_Doctors_Visited),
    Median_kunjungan      = median(Num_Doctors_Visited),
    Variansi_kunjungan    = var(Num_Doctors_Visited),
    Min_kunjungan         = min(Num_Doctors_Visited),
    Max_kunjungan         = max(Num_Doctors_Visited)
  )
knitr::kable(ringkasan_poisson, caption = "Ringkasan Data Poisson (Jumlah Kunjungan Dokter)",
             digits = 3)
Ringkasan Data Poisson (Jumlah Kunjungan Dokter)
Jumlah_observasi Rata_rata_kunjungan Median_kunjungan Variansi_kunjungan Min_kunjungan Max_kunjungan
714 2.112 2 0.467 1 3
dist_poisson <- data_poisson %>%
  dplyr::count(Num_Doctors_Visited) %>%
  dplyr::mutate(Persen = scales::percent(n / sum(n)))
knitr::kable(dist_poisson, caption = "Distribusi Jumlah Dokter yang Dikunjungi")
Distribusi Jumlah Dokter yang Dikunjungi
Num_Doctors_Visited n Persen
1 131 18%
2 372 52%
3 211 30%
ggplot(data_poisson, aes(x = Num_Doctors_Visited)) +
  geom_histogram(binwidth = 1, fill = "#2f7f73", color = "white") +
  labs(
    title = "Distribusi Jumlah Dokter yang Dikunjungi Pasien Lansia",
    x     = "Jumlah dokter dikunjungi",
    y     = "Jumlah pasien"
  ) +
  theme_minimal()

Jika rata-rata dan variansi jumlah kunjungan relatif dekat, maka asumsi Poisson (mean = variansi) terpenuhi dengan baik. Jika variansi jauh lebih besar dari rata-rata, perlu dipertimbangkan quasi-Poisson atau Negative Binomial sebagai alternatif.

9.5 Pembentukan Model Regresi Poisson

Model Poisson yang digunakan untuk data individual tanpa offset:

\[ \log(\mu_i) = \beta_0 + \beta_1\,\text{Age}_i + \beta_2\,\text{Physical\_Health}_i + \cdots + \beta_p X_{pi} \]

Karena data berbentuk individual (satu baris per pasien), model Poisson tanpa offset digunakan untuk mengestimasi rata-rata jumlah dokter yang dikunjungi (\(\mu_i\)) secara langsung.

model_poisson <- glm(
  Num_Doctors_Visited ~ Age + Physical_Health + Mental_Health + Dental_Health +
    Employment + Stress_Sleep + Medication_Sleep + Pain_Sleep +
    Bathroom_Sleep + Unknown_Sleep + Trouble_Sleeping +
    Prescription_Sleep_Med + Race + Gender,
  data   = data_poisson,
  family = poisson(link = "log")
)

summary(model_poisson)
## 
## Call:
## glm(formula = Num_Doctors_Visited ~ Age + Physical_Health + Mental_Health + 
##     Dental_Health + Employment + Stress_Sleep + Medication_Sleep + 
##     Pain_Sleep + Bathroom_Sleep + Unknown_Sleep + Trouble_Sleeping + 
##     Prescription_Sleep_Med + Race + Gender, family = poisson(link = "log"), 
##     data = data_poisson)
## 
## Coefficients: (1 not defined because of singularities)
##                          Estimate Std. Error z value Pr(>|z|)  
## (Intercept)             0.6157139  0.2493428   2.469   0.0135 *
## Age                            NA         NA      NA       NA  
## Physical_Health         0.0644736  0.0348218   1.852   0.0641 .
## Mental_Health          -0.0148086  0.0305598  -0.485   0.6280  
## Dental_Health          -0.0126340  0.0209388  -0.603   0.5463  
## Employment              0.0362067  0.0465152   0.778   0.4363  
## Stress_Sleep            0.0515882  0.0652498   0.791   0.4292  
## Medication_Sleep        0.1333316  0.1086001   1.228   0.2195  
## Pain_Sleep              0.0061158  0.0679182   0.090   0.9282  
## Bathroom_Sleep          0.0416547  0.0566048   0.736   0.4618  
## Unknown_Sleep           0.0290452  0.0614102   0.473   0.6362  
## Trouble_Sleeping        0.0130477  0.0422874   0.309   0.7577  
## Prescription_Sleep_Med -0.0511115  0.0468408  -1.091   0.2752  
## Race                   -0.0200849  0.0264055  -0.761   0.4469  
## Gender                  0.0009173  0.0528779   0.017   0.9862  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 167.31  on 713  degrees of freedom
## Residual deviance: 157.05  on 700  degrees of freedom
## AIC: 2050.6
## 
## Number of Fisher Scoring iterations: 4

9.6 Persamaan Model Hasil Estimasi Poisson

Chunk model-poisson di atas sudah membangun model regresi Poisson menggunakan glm() dengan family = poisson(link = "log"). Bentuk hasil estimasinya adalah:

\[ \log(\hat{\mu_i}) = \hat{\beta}_0 + \hat{\beta}_1 X_{1i} + \cdots + \hat{\beta}_p X_{pi} \]

Koefisien yang dieksponensialkan menjadi Incidence Rate Ratio (IRR). IRR > 1 berarti variabel tersebut berkaitan dengan rata-rata jumlah kunjungan dokter yang lebih tinggi.

Parameter yang Membentuk Persamaan Regresi Poisson
Variabel_Mudah term estimate IRR
(Intercept) (Intercept) 0.6157 1.8510
kelompok usia pasien Age NA NA
Physical_Health Physical_Health 0.0645 1.0666
Mental_Health Mental_Health -0.0148 0.9853
Dental_Health Dental_Health -0.0126 0.9874
status pekerjaan pasien Employment 0.0362 1.0369
Stress_Sleep Stress_Sleep 0.0516 1.0529
Medication_Sleep Medication_Sleep 0.1333 1.1426
Pain_Sleep Pain_Sleep 0.0061 1.0061
Bathroom_Sleep Bathroom_Sleep 0.0417 1.0425
Unknown_Sleep Unknown_Sleep 0.0290 1.0295
Trouble_Sleeping Trouble_Sleeping 0.0130 1.0131
Prescription_Sleep_Med Prescription_Sleep_Med -0.0511 0.9502
kategori ras pasien Race -0.0201 0.9801
jenis kelamin pasien Gender 0.0009 1.0009

9.7 Koefisien dan Incidence Rate Ratio Model Poisson

tabel_poisson <- broom::tidy(model_poisson) %>%
  dplyr::mutate(
    Incidence_Rate_Ratio = exp(estimate),
    Perubahan_Persen     = 100 * (Incidence_Rate_Ratio - 1),
    Variabel_Mudah       = purrr::map_chr(term, label_poisson),
    Keputusan_5persen    = ifelse(p.value < 0.05, "Signifikan", "Tidak signifikan")
  ) %>%
  dplyr::arrange(p.value)

knitr::kable(
  tabel_poisson %>%
    dplyr::select(Variabel_Mudah, term, estimate, std.error, statistic,
                  p.value, Incidence_Rate_Ratio, Perubahan_Persen, Keputusan_5persen),
  digits  = 4,
  caption = "Koefisien dan Incidence Rate Ratio Regresi Poisson"
)
Koefisien dan Incidence Rate Ratio Regresi Poisson
Variabel_Mudah term estimate std.error statistic p.value Incidence_Rate_Ratio Perubahan_Persen Keputusan_5persen
(Intercept) (Intercept) 0.6157 0.2493 2.4693 0.0135 1.8510 85.0978 Signifikan
Physical_Health Physical_Health 0.0645 0.0348 1.8515 0.0641 1.0666 6.6597 Tidak signifikan
Medication_Sleep Medication_Sleep 0.1333 0.1086 1.2277 0.2195 1.1426 14.2629 Tidak signifikan
Prescription_Sleep_Med Prescription_Sleep_Med -0.0511 0.0468 -1.0912 0.2752 0.9502 -4.9827 Tidak signifikan
Stress_Sleep Stress_Sleep 0.0516 0.0652 0.7906 0.4292 1.0529 5.2942 Tidak signifikan
status pekerjaan pasien Employment 0.0362 0.0465 0.7784 0.4363 1.0369 3.6870 Tidak signifikan
kategori ras pasien Race -0.0201 0.0264 -0.7606 0.4469 0.9801 -1.9885 Tidak signifikan
Bathroom_Sleep Bathroom_Sleep 0.0417 0.0566 0.7359 0.4618 1.0425 4.2534 Tidak signifikan
Dental_Health Dental_Health -0.0126 0.0209 -0.6034 0.5463 0.9874 -1.2554 Tidak signifikan
Mental_Health Mental_Health -0.0148 0.0306 -0.4846 0.6280 0.9853 -1.4699 Tidak signifikan
Unknown_Sleep Unknown_Sleep 0.0290 0.0614 0.4730 0.6362 1.0295 2.9471 Tidak signifikan
Trouble_Sleeping Trouble_Sleeping 0.0130 0.0423 0.3085 0.7577 1.0131 1.3133 Tidak signifikan
Pain_Sleep Pain_Sleep 0.0061 0.0679 0.0900 0.9282 1.0061 0.6135 Tidak signifikan
jenis kelamin pasien Gender 0.0009 0.0529 0.0173 0.9862 1.0009 0.0918 Tidak signifikan
kelompok usia pasien Age NA NA NA NA NA NA NA

9.8 Evaluasi Model Poisson

AIC_poisson <- AIC(model_poisson)
AIC_poisson
## [1] 2050.627
dispersion_poisson <- deviance(model_poisson) / df.residual(model_poisson)
dispersion_poisson
## [1] 0.2243558
if (dispersion_poisson > 1.5) {
  cat("Terdapat indikasi overdispersion (rasio deviance/df =",
      round(dispersion_poisson, 3),
      "). Inferensi p-value sebaiknya dilihat juga menggunakan quasi-Poisson.\n")
} else {
  cat("Tidak ada indikasi overdispersion yang kuat berdasarkan rasio deviance/df residual =",
      round(dispersion_poisson, 3), ".\n")
}
## Tidak ada indikasi overdispersion yang kuat berdasarkan rasio deviance/df residual = 0.224 .

9.9 Pemeriksaan Overdispersion dengan Pearson Dispersion

Selain rasio deviance terhadap derajat bebas residual, overdispersion juga dapat diperiksa menggunakan Pearson dispersion.

\[ \hat{\phi} = \frac{\sum r_{P,i}^2}{df} \]

Jika nilai dispersion jauh lebih besar dari 1, maka terdapat indikasi overdispersion.

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

tibble::tibble(
  Ukuran        = "Pearson dispersion",
  Nilai         = pearson_dispersion,
  Interpretasi  = dplyr::case_when(
    pearson_dispersion < 1.5 ~ "Tidak ada indikasi overdispersion berat",
    pearson_dispersion < 2.5 ~ "Ada indikasi overdispersion sedang",
    TRUE                     ~ "Ada indikasi overdispersion kuat"
  )
) %>%
  dplyr::mutate(Nilai = round(Nilai, 3)) %>%
  knitr::kable(caption = "Pemeriksaan Overdispersion dengan Pearson Dispersion")
Pemeriksaan Overdispersion dengan Pearson Dispersion
Ukuran Nilai Interpretasi
Pearson dispersion 0.212 Tidak ada indikasi overdispersion berat

Jika nilai Pearson dispersion jauh lebih besar dari 1, maka variasi jumlah kunjungan dokter lebih besar daripada yang diasumsikan oleh model Poisson. Dalam kondisi tersebut, hasil signifikansi dari model Poisson standar perlu dibaca hati-hati dan model quasi-Poisson lebih aman digunakan sebagai pembanding.

model_quasipoisson <- glm(
  Num_Doctors_Visited ~ Age + Physical_Health + Mental_Health + Dental_Health +
    Employment + Stress_Sleep + Medication_Sleep + Pain_Sleep +
    Bathroom_Sleep + Unknown_Sleep + Trouble_Sleeping +
    Prescription_Sleep_Med + Race + Gender,
  data   = data_poisson,
  family = quasipoisson(link = "log")
)

tabel_quasipoisson <- broom::tidy(model_quasipoisson) %>%
  dplyr::mutate(
    Incidence_Rate_Ratio = exp(estimate),
    Perubahan_Persen     = 100 * (Incidence_Rate_Ratio - 1),
    Variabel_Mudah       = purrr::map_chr(term, label_poisson),
    Keputusan_5persen    = ifelse(p.value < 0.05, "Signifikan", "Tidak signifikan")
  ) %>%
  dplyr::arrange(p.value)

knitr::kable(
  tabel_quasipoisson %>%
    dplyr::select(Variabel_Mudah, term, estimate, std.error, statistic,
                  p.value, Incidence_Rate_Ratio, Perubahan_Persen, Keputusan_5persen),
  digits  = 4,
  caption = "Koefisien dan IRR Quasi-Poisson"
)
Koefisien dan IRR Quasi-Poisson
Variabel_Mudah term estimate std.error statistic p.value Incidence_Rate_Ratio Perubahan_Persen Keputusan_5persen
(Intercept) (Intercept) 0.6157 0.1149 5.3599 0.0000 1.8510 85.0978 Signifikan
Physical_Health Physical_Health 0.0645 0.0160 4.0189 0.0001 1.0666 6.6597 Signifikan
Medication_Sleep Medication_Sleep 0.1333 0.0500 2.6649 0.0079 1.1426 14.2629 Signifikan
Prescription_Sleep_Med Prescription_Sleep_Med -0.0511 0.0216 -2.3685 0.0181 0.9502 -4.9827 Signifikan
Stress_Sleep Stress_Sleep 0.0516 0.0301 1.7161 0.0866 1.0529 5.2942 Tidak signifikan
status pekerjaan pasien Employment 0.0362 0.0214 1.6895 0.0916 1.0369 3.6870 Tidak signifikan
kategori ras pasien Race -0.0201 0.0122 -1.6510 0.0992 0.9801 -1.9885 Tidak signifikan
Bathroom_Sleep Bathroom_Sleep 0.0417 0.0261 1.5973 0.1106 1.0425 4.2534 Tidak signifikan
Dental_Health Dental_Health -0.0126 0.0096 -1.3097 0.1907 0.9874 -1.2554 Tidak signifikan
Mental_Health Mental_Health -0.0148 0.0141 -1.0518 0.2932 0.9853 -1.4699 Tidak signifikan
Unknown_Sleep Unknown_Sleep 0.0290 0.0283 1.0266 0.3050 1.0295 2.9471 Tidak signifikan
Trouble_Sleeping Trouble_Sleeping 0.0130 0.0195 0.6697 0.5033 1.0131 1.3133 Tidak signifikan
Pain_Sleep Pain_Sleep 0.0061 0.0313 0.1955 0.8451 1.0061 0.6135 Tidak signifikan
jenis kelamin pasien Gender 0.0009 0.0244 0.0377 0.9700 1.0009 0.0918 Tidak signifikan
kelompok usia pasien Age NA NA NA NA NA NA NA
data_poisson$fitted_poisson <- fitted(model_poisson)

ggplot(data_poisson, aes(x = fitted_poisson, y = Num_Doctors_Visited)) +
  geom_point(alpha = 0.35, color = "#2f7f73") +
  geom_abline(slope = 1, intercept = 0, linetype = "dashed", color = "#b84f5a") +
  labs(
    title = "Observed vs Fitted Jumlah Dokter yang Dikunjungi",
    x     = "Nilai fitted model",
    y     = "Jumlah kunjungan dokter aktual"
  ) +
  theme_minimal()

9.10 Visualisasi Prediksi Rata-rata Kunjungan Dokter

Visualisasi ini digunakan untuk melihat prediksi rata-rata kunjungan dokter berdasarkan status kesehatan fisik dan mental pasien.

grid_poisson <- expand.grid(
  Physical_Health       = 1:5,
  Mental_Health         = c(1, 3, 5),
  Age                   = median(data_poisson$Age, na.rm = TRUE),
  Dental_Health         = median(data_poisson$Dental_Health, na.rm = TRUE),
  Employment            = modus(data_poisson$Employment),
  Stress_Sleep          = modus(data_poisson$Stress_Sleep),
  Medication_Sleep      = modus(data_poisson$Medication_Sleep),
  Pain_Sleep            = modus(data_poisson$Pain_Sleep),
  Bathroom_Sleep        = modus(data_poisson$Bathroom_Sleep),
  Unknown_Sleep         = modus(data_poisson$Unknown_Sleep),
  Trouble_Sleeping      = modus(data_poisson$Trouble_Sleeping),
  Prescription_Sleep_Med= modus(data_poisson$Prescription_Sleep_Med),
  Race                  = modus(data_poisson$Race),
  Gender                = modus(data_poisson$Gender)
)

pred_rate_poisson <- predict(
  model_poisson,
  newdata = grid_poisson,
  type    = "link",
  se.fit  = TRUE
)

grid_poisson_plot <- grid_poisson %>%
  dplyr::mutate(
    fit_link        = pred_rate_poisson$fit,
    se_link         = pred_rate_poisson$se.fit,
    Predicted_Count = exp(fit_link),
    Count_Lower     = exp(fit_link - 1.96 * se_link),
    Count_Upper     = exp(fit_link + 1.96 * se_link),
    Mental_Health   = factor(Mental_Health,
                             labels = c("Mental baik (1)", "Mental sedang (3)", "Mental buruk (5)"))
  )

ggplot(
  grid_poisson_plot,
  aes(x = Physical_Health, y = Predicted_Count, color = Mental_Health, group = Mental_Health)
) +
  geom_line(linewidth = 1.1) +
  geom_ribbon(aes(ymin = Count_Lower, ymax = Count_Upper, fill = Mental_Health),
              alpha = 0.15, color = NA) +
  labs(
    title    = "Prediksi Rata-rata Jumlah Kunjungan Dokter",
    subtitle = "Berdasarkan status kesehatan fisik dan mental pasien",
    x        = "Status kesehatan fisik (1 = sangat baik, 5 = sangat buruk)",
    y        = "Prediksi rata-rata jumlah kunjungan dokter",
    color    = "Status mental",
    fill     = "Status mental"
  ) +
  theme_minimal()

Grafik ini menunjukkan prediksi rata-rata jumlah kunjungan dokter menurut status kesehatan fisik dan mental pasien, dengan variabel lain ditahan pada nilai paling umum. Jika garis cenderung naik ketika status kesehatan fisik memburuk (nilai lebih tinggi), maka model menunjukkan bahwa kondisi fisik yang lebih buruk berkaitan dengan lebih banyak kunjungan dokter.

9.11 Interpretasi dan Kesimpulan Model Poisson

if (dispersion_poisson > 1.5) {
  tabel_pakai_poisson <- tabel_quasipoisson
  metode_pakai        <- "quasi-Poisson"
} else {
  tabel_pakai_poisson <- tabel_poisson
  metode_pakai        <- "Poisson"
}

sig_poisson <- tabel_pakai_poisson %>%
  dplyr::filter(term != "(Intercept)", p.value < 0.05) %>%
  dplyr::arrange(p.value)

cat("Rasio deviance/df residual adalah **", fmt_num(dispersion_poisson, 3), "**. ", sep = "")

Rasio deviance/df residual adalah 0.224.

if (dispersion_poisson > 1.5) {
  cat("Nilai ini menunjukkan overdispersion, sehingga interpretasi signifikansi lebih aman mengacu pada model quasi-Poisson.\n\n")
} else {
  cat("Nilai ini tidak menunjukkan overdispersion yang kuat, sehingga model Poisson standar masih memadai.\n\n")
}

Nilai ini tidak menunjukkan overdispersion yang kuat, sehingga model Poisson standar masih memadai.

cat("Dengan acuan inferensi ", metode_pakai, ", terdapat ", nrow(sig_poisson),
    " prediktor signifikan pada taraf 5%.\n\n", sep = "")

Dengan acuan inferensi Poisson, terdapat 0 prediktor signifikan pada taraf 5%.

if (nrow(sig_poisson) > 0) {
  cat("Prediktor paling penting berdasarkan p-value terkecil adalah:\n\n")
  top_poisson <- head(sig_poisson, 10)
  for (i in seq_len(nrow(top_poisson))) {
    row  <- top_poisson[i, ]
    arah <- ifelse(row$Incidence_Rate_Ratio > 1,
                   "berkaitan dengan rata-rata kunjungan dokter yang lebih tinggi",
                   "berkaitan dengan rata-rata kunjungan dokter yang lebih rendah")
    cat("- **", row$Variabel_Mudah, "** memiliki IRR = ", fmt_num(row$Incidence_Rate_Ratio, 3),
        " dengan perubahan sekitar ", fmt_num(row$Perubahan_Persen, 2), "%",
        " dan p-value = ", fmt_p(row$p.value), ", sehingga variabel ini ", arah,
        ", dengan asumsi variabel lain tetap.\n", sep = "")
  }
  cat("\n")
}

cat("AIC model Poisson adalah **", fmt_num(AIC_poisson, 2), "**.\n\n", sep = "")

AIC model Poisson adalah 2,050.63.

cat("**Kesimpulan substantif:** jumlah dokter yang dikunjungi pasien lansia tidak hanya dipengaruhi oleh status kesehatan fisik dan mental, tetapi juga oleh gangguan tidur, status pekerjaan, dan karakteristik demografi. Karena data berbentuk individual tanpa offset, model Poisson mengestimasi rata-rata jumlah kunjungan dokter secara langsung berdasarkan kondisi tiap pasien.\n")

Kesimpulan substantif: jumlah dokter yang dikunjungi pasien lansia tidak hanya dipengaruhi oleh status kesehatan fisik dan mental, tetapi juga oleh gangguan tidur, status pekerjaan, dan karakteristik demografi. Karena data berbentuk individual tanpa offset, model Poisson mengestimasi rata-rata jumlah kunjungan dokter secara langsung berdasarkan kondisi tiap pasien.

10 Perbandingan Empat Model

Model Dataset Respon Fungsi link Ukuran interpretasi Evaluasi utama
Logistik biner Telco Customer Churn Churn 0/1 Logit Odds ratio Confusion matrix, AIC, ROC/AUC
Ordinal Student Performance (Portugis) Performa G3 (Rendah/Sedang/Tinggi) Cumulative logit MASS::polr: logit(P(Y <= j)) = alpha_j - x beta exp(beta) sebagai OR menuju performa lebih tinggi Confusion matrix, AIC, proportional odds check
Multinomial HR Employee Attrition JobRole (9 kategori) Baseline-category logit vs Sales Executive Relative risk ratio Confusion matrix, AIC, prediksi probabilitas, catatan IIA
Poisson NPHA Doctor Visits Jumlah kunjungan dokter Log (tanpa offset, data individu) Incidence rate ratio dan perubahan persen AIC, deviance dispersion, Pearson dispersion, overdispersion

11 Kesimpulan Umum

cat("Berdasarkan empat analisis yang dilakukan, setiap model menjawab pertanyaan substantif yang berbeda dari empat konteks data yang beragam.\n\n")

Berdasarkan empat analisis yang dilakukan, setiap model menjawab pertanyaan substantif yang berbeda dari empat konteks data yang beragam.

cat("- **Model biner** menjelaskan peluang pelanggan telekomunikasi melakukan churn (berhenti berlangganan). Model ini cocok karena respon hanya terdiri dari dua kategori: churn dan tidak churn. Faktor seperti jenis kontrak, lama berlangganan, dan layanan internet menjadi penentu utama.\n")
  • Model biner menjelaskan peluang pelanggan telekomunikasi melakukan churn (berhenti berlangganan). Model ini cocok karena respon hanya terdiri dari dua kategori: churn dan tidak churn. Faktor seperti jenis kontrak, lama berlangganan, dan layanan internet menjadi penentu utama.
cat("- **Model ordinal** menjelaskan tingkat performa akademik siswa berdasarkan nilai akhir G3 yang dikelompokkan menjadi Rendah, Sedang, dan Tinggi. Model ini cocok karena performa memiliki urutan yang bermakna. Faktor seperti kegagalan kelas sebelumnya, waktu belajar, dan konsumsi alkohol menjadi prediktor penting.\n")
  • Model ordinal menjelaskan tingkat performa akademik siswa berdasarkan nilai akhir G3 yang dikelompokkan menjadi Rendah, Sedang, dan Tinggi. Model ini cocok karena performa memiliki urutan yang bermakna. Faktor seperti kegagalan kelas sebelumnya, waktu belajar, dan konsumsi alkohol menjadi prediktor penting.
cat("- **Model multinomial** menjelaskan peran jabatan karyawan (JobRole) di perusahaan. Model ini cocok karena jabatan memiliki beberapa kategori yang tidak diperlakukan sebagai urutan. Pendapatan bulanan, departemen, dan total pengalaman kerja membedakan antar peran jabatan.\n")
  • Model multinomial menjelaskan peran jabatan karyawan (JobRole) di perusahaan. Model ini cocok karena jabatan memiliki beberapa kategori yang tidak diperlakukan sebagai urutan. Pendapatan bulanan, departemen, dan total pengalaman kerja membedakan antar peran jabatan.
cat("- **Model Poisson** menjelaskan jumlah dokter yang dikunjungi oleh pasien lansia. Model ini cocok karena respon berupa hitungan kejadian yang bersifat diskret dan non-negatif. Status kesehatan fisik, mental, dan gangguan tidur menjadi faktor utama penjelas.\n\n")
  • Model Poisson menjelaskan jumlah dokter yang dikunjungi oleh pasien lansia. Model ini cocok karena respon berupa hitungan kejadian yang bersifat diskret dan non-negatif. Status kesehatan fisik, mental, dan gangguan tidur menjadi faktor utama penjelas.
cat("Secara keseluruhan, analisis ini menunjukkan bahwa pemilihan model harus disesuaikan dengan skala pengukuran variabel respon. Model biner untuk respon dikotomi, model ordinal untuk respon bertingkat, model multinomial untuk respon nominal banyak kategori, dan model Poisson untuk data cacah. Keempat model ini merupakan bagian penting dari keluarga generalized linear model (GLM).\n")

Secara keseluruhan, analisis ini menunjukkan bahwa pemilihan model harus disesuaikan dengan skala pengukuran variabel respon. Model biner untuk respon dikotomi, model ordinal untuk respon bertingkat, model multinomial untuk respon nominal banyak kategori, dan model Poisson untuk data cacah. Keempat model ini merupakan bagian penting dari keluarga generalized linear model (GLM).

12 Sumber Data

  1. https://www.kaggle.com/datasets/blastchar/telco-customer-churn (Biner — Telco Customer Churn)
  2. https://archive.ics.uci.edu/dataset/320/student+performance (Ordinal — Student Performance Portugis)
  3. https://www.kaggle.com/datasets/pavansubhasht/ibm-hr-analytics-attrition-dataset (Multinomial — HR Employee Attrition)
  4. https://www.icpsr.umich.edu/web/NACDA/studies/37305 (Poisson — NPHA Doctor Visits)