Klasifikasi kelulusan seleksi masuk perguruan tinggi dengan k-Nearest Neighbors

# Contoh kode R: k-Nearest Neighbors
suppressPackageStartupMessages({
  library(ggplot2)
  library(patchwork)
  library(class)
  library(pROC)
})
set.seed(42)
# ---- Gaya grafik ----
pal <- c("#1f6feb", "#e8590c", "#2f9e44", "#ae3ec9", "#f59f00")
warna <- list(bg = "#f4f6fa", panel = "#ffffff", teks = "#1f2933", redup = "#5f6b78", grid = "#e3e8ee", bar = "#1d3557")
tata <- 1
theme_set(theme_minimal(base_size = 11) + theme(legend.position = "bottom"))

# ---- Data sintetis: kelas biner dari dua fitur utama dan satu fitur tambahan ----
n <- 200
positif <- function(x, m) if (m > 0) pmax(x, 0.05 * m) else x   # cegah nilai negatif yang tidak wajar
d <- data.frame(academic_aptitude_score = positif(rnorm(n, 50, 10), 50),
                mean_report_card_grade = positif(rnorm(n, 30, 6), 30),
                x3 = round(rnorm(n, 50, 10), 2))
z1 <- as.numeric(scale(d$academic_aptitude_score))
z2 <- as.numeric(scale(d$mean_report_card_grade))
z3 <- as.numeric(scale(d$x3))
arah <- sample(c(-1, 1), 2, replace = TRUE)
skor <- arah[1] * runif(1, 0.8, 1.3) * z1 + arah[2] * runif(1, 0.6, 1.1) * z2 +
  sample(c(-1, 1), 1) * runif(1, 0.5, 0.9) * z1 * z2 + runif(1, 0.1, 0.35) * z3
skor <- skor - quantile(skor, runif(1, 0.4, 0.6))
d$admission_outcome <- rbinom(n, 1, plogis(runif(1, 2.5, 4) * skor))
write.csv(d, "data_contoh.csv", row.names = FALSE)

# ---- Pembagian data latih/uji (70:30) dan standarisasi memakai rerata/SB data latih ----
d$kelas <- factor(d$admission_outcome, levels = c(0, 1), labels = c("No", "Yes"))
fitur <- c("academic_aptitude_score", "mean_report_card_grade", "x3")
idx <- unlist(lapply(split(seq_len(n), d$kelas), function(i) sample(i, round(0.7 * length(i)))))
latih <- d[idx, ]
uji <- d[-idx, ]
rerata <- sapply(latih[fitur], mean)
sb <- sapply(latih[fitur], sd)
xl <- scale(latih[fitur], center = rerata, scale = sb)
xu <- scale(uji[fitur], center = rerata, scale = sb)

# ---- Pemilihan k: validasi silang leave-one-out (class::knn.cv) untuk k ganjil 1..25 ----
daftar_k <- seq(1, 25, by = 2)
akurasi_cv <- sapply(daftar_k, function(k) mean(class::knn.cv(xl, latih$kelas, k = k) == latih$kelas))
akurasi_uji_k <- sapply(daftar_k, function(k) mean(class::knn(xl, xu, latih$kelas, k = k) == uji$kelas))
k_best <- daftar_k[which.max(akurasi_cv)]
print(data.frame(k = daftar_k, akurasi_cv = round(akurasi_cv, 3)))
cat("k terbaik (LOOCV):", k_best, "\n")

# ---- Model akhir dan evaluasi pada data uji ----
pr <- class::knn(xl, xu, latih$kelas, k = k_best, prob = TRUE)
p_menang <- attr(pr, "prob")                       # proporsi suara kelas pemenang
prob <- ifelse(pr == "Yes", p_menang, 1 - p_menang) # peluang kelas Yes
konf <- table(Actual = uji$kelas, Predicted = pr)
print(konf)
tp <- konf["Yes", "Yes"]; fp <- konf["No", "Yes"]; fn <- konf["Yes", "No"]; tn <- konf["No", "No"]
akurasi <- sum(diag(konf)) / sum(konf)
presisi <- if (tp + fp > 0) tp / (tp + fp) else 0
recall <- if (tp + fn > 0) tp / (tp + fn) else 0
spesifisitas <- if (tn + fp > 0) tn / (tn + fp) else 0
f1 <- if (presisi + recall > 0) 2 * presisi * recall / (presisi + recall) else 0
roc_knn <- pROC::roc(uji$kelas, prob, levels = c("No", "Yes"), direction = "<", quiet = TRUE)
auc_knn <- as.numeric(pROC::auc(roc_knn))
cat(sprintf("Akurasi = %.3f, F1 = %.3f, AUC = %.3f\n", akurasi, f1, auc_knn))


# ---- Grafik ----
dk <- rbind(data.frame(k = daftar_k, akurasi = akurasi_cv, jenis = "LOOCV (training)"),
            data.frame(k = daftar_k, akurasi = akurasi_uji_k, jenis = "Test set"))
p1 <- ggplot(dk, aes(k, akurasi, colour = jenis)) +
  geom_vline(xintercept = k_best, linetype = "dashed", colour = warna$redup) +
  geom_line(linewidth = 0.9) +
  geom_point(size = 1.8) +
  annotate("text", x = k_best, y = min(dk$akurasi), label = paste0(" best k = ", k_best), hjust = 0, vjust = 0,
           colour = warna$teks, size = 3.1) +
  scale_colour_manual(values = c(pal[1], pal[2]), name = NULL) +
  scale_x_continuous(breaks = seq(1, 25, by = 4)) +
  labs(title = "Accuracy vs number of neighbours", x = "k", y = "Accuracy") +
  theme(legend.position = "bottom")
gx <- seq(min(d$academic_aptitude_score), max(d$academic_aptitude_score), length.out = 90)
gy <- seq(min(d$mean_report_card_grade), max(d$mean_report_card_grade), length.out = 90)
bidang <- expand.grid(academic_aptitude_score = gx, mean_report_card_grade = gy)
bidang$x3 <- rerata["x3"]                           # fitur tambahan ditahan pada reratanya
xb <- scale(bidang[fitur], center = rerata, scale = sb)
pb <- class::knn(xl, xb, latih$kelas, k = k_best, prob = TRUE)
bidang$prob <- ifelse(pb == "Yes", attr(pb, "prob"), 1 - attr(pb, "prob"))
p2 <- ggplot() +
  geom_raster(data = bidang, aes(academic_aptitude_score, mean_report_card_grade, fill = prob), alpha = 0.5) +
  geom_contour(data = bidang, aes(academic_aptitude_score, mean_report_card_grade, z = prob), breaks = 0.5, colour = warna$teks, linewidth = 0.6) +
  geom_point(data = uji, aes(academic_aptitude_score, mean_report_card_grade, colour = kelas, shape = kelas), size = 1.9) +
  scale_fill_gradient2(low = pal[1], mid = warna$panel, high = pal[2], midpoint = 0.5, guide = "none") +
  scale_colour_manual(values = c(No = pal[1], Yes = pal[2]), name = "admission_outcome") +
  scale_shape_manual(values = c(No = 16, Yes = 17), name = "admission_outcome") +
  coord_cartesian(expand = FALSE) +
  labs(title = sprintf("Decision regions, k = %d (test points)", k_best), subtitle = "Shading = share of Yes neighbours",
       x = "academic_aptitude_score", y = "mean_report_card_grade")
dkf <- as.data.frame(konf)
dkf$persen <- dkf$Freq / sum(dkf$Freq)
dkf$benar <- dkf$Actual == dkf$Predicted
p3 <- ggplot(dkf, aes(Predicted, Actual)) +
  geom_tile(aes(fill = benar, alpha = persen), colour = warna$panel, linewidth = 1) +
  geom_text(aes(label = sprintf("%d\n(%.0f%%)", Freq, 100 * persen)), colour = warna$teks, size = 3.4, lineheight = 0.9) +
  scale_fill_manual(values = c(`TRUE` = pal[3], `FALSE` = pal[5]), guide = "none") +
  scale_alpha(range = c(0.25, 0.85), guide = "none") +
  scale_y_discrete(limits = c("Yes", "No")) +
  labs(title = "Confusion matrix (test data)", x = "Predicted", y = "Actual") +
  theme(panel.grid = element_blank())
p_hasil <- if (tata == 1) (p2 | (p1 / p3)) + plot_layout(widths = c(1.35, 1)) else (p2 / (p1 | p3)) + plot_layout(heights = c(1.3, 1))

print(p_hasil)