Penggunaan Metrik Evaluasi AIC-CV dalam Menentukan Model Terbaik: Kasus CCPP
Dosen Pembimbing:
Prof. Dr. I Gede Nyoman Mindra Jaya, M.Si.
Disusun Oleh:
Rasendriya Nandana Kurniawan - 140720260011
Yuda Taufiqurahman Wenske - 140720260017
PROGRAM STUDI MAGISTER STATISTIKA TERAPAN
DEPARTEMEN STATISTIKA
FAKULTAS MATEMATIKA DAN ILMU PENGETAHUAN ALAM
UNIVERSITAS PADJADJARAN
TAHUN AJARAN 2025/2026
Bab 1. Pendahuluan
1.1. Latar Belakang
Energi listrik merupakan kebutuhan dasar yang terus meningkat seiring pertumbuhan penduduk, industrialisasi, dan digitalisasi. Dalam sistem pembangkit listrik, khususnya Combined Cycle Power Plant (CCPP), kemampuan memprediksi daya keluaran listrik (power output) secara akurat menjadi faktor penting untuk meningkatkan efisiensi operasional, menekan biaya bahan bakar, dan menjaga keandalan pasokan. CCPP bekerja dengan menggabungkan turbin gas dan turbin uap, sehingga kinerjanya dipengaruhi oleh berbagai kondisi lingkungan seperti suhu udara, kelembapan, tekanan udara, dan vakum buang.
Dataset yang digunakan dalam penelitian ini memuat lima variabel utama, yaitu AT (suhu lingkungan), V (vakum buang), AP (tekanan udara ambien), RH (kelembapan relatif), dan PE (daya keluaran listrik). Variabel AT, V, AP, dan RH berperan sebagai prediktor, sedangkan PE merupakan variabel respons yang akan diprediksi. Hubungan antara prediktor dan respons tidak selalu linear, sehingga pendekatan regresi linear sederhana berpotensi kurang mampu menangkap pola yang kompleks.
Perkembangan metode statistika dan pembelajaran mesin memberikan alternatif untuk memodelkan hubungan non-linear dan multikolinearitas. Polynomial Regression dapat menangkap efek non-linear melalui suku pangkat, Elastic Net mampu menangani multikolinearitas dan melakukan seleksi variabel melalui regularisasi, sedangkan Generalized Additive Model (GAM) memodelkan hubungan non-linear secara fleksibel melalui fungsi halus (smooth functions). Ketiga metode tersebut memiliki karakteristik yang berbeda, sehingga perlu dibandingkan secara objektif.
Pemilihan model yang tepat memerlukan kriteria yang dapat menyeimbangkan kecocokan model dan kompleksitas. AICc (Corrected Akaike Information Criterion) digunakan untuk pemilihan model pada data berukuran kecil hingga sedang, sedangkan 10-fold Cross-Validation digunakan untuk mengukur kemampuan generalisasi model terhadap data yang tidak digunakan dalam pelatihan. Selain itu, dataset dalam penelitian ini tersedia dalam beberapa sheet, sehingga perlu dianalisis apakah kinerja model berbeda secara signifikan antar sheet.
Berdasarkan uraian tersebut, penelitian ini dilakukan untuk membandingkan Polynomial Regression, Elastic Net, dan GAM dalam memprediksi PE pada data CCPP, memilih model terbaik berdasarkan AICc dan 10-fold Cross-Validation, serta menganalisis perbedaan kinerja model antar sheet menggunakan ANOVA dan Kruskal-Wallis. Hasil penelitian diharapkan dapat memberikan gambaran mengenai metode yang paling sesuai dan memberikan rekomendasi bagi pengelolaan data energi.
1.2. Identifikasi Masalah
Berdasarkan latar belakang, masalah dalam penelitian ini dapat diidentifikasi sebagai berikut.
- Prediksi PE pada CCPP penting, tetapi hubungan antara AT, V, AP, RH dan PE tidak selalu linear.
- Terdapat potensi multikolinearitas dan interaksi antar prediktor yang dapat memengaruhi kinerja model regresi.
- Metode Polynomial Regression, Elastic Net, dan GAM memiliki kelebihan dan keterbatasan yang berbeda, sehingga perlu dibandingkan secara empiris.
- Pemilihan model terbaik tidak cukup hanya berdasarkan satu kriteria, sehingga diperlukan perbandingan AICc dan 10-fold Cross-Validation.
- Dataset terdiri atas beberapa sheet, sehingga perlu diketahui apakah kinerja model berbeda nyata antar sheet.
- Belum banyak penelitian yang secara khusus membandingkan ketiga metode tersebut pada dataset CCPP dengan pendekatan pemilihan model AICc dan cross-validation.
1.3. Tujuan
Penelitian ini dimaksudkan untuk menerapkan dan membandingkan metode Polynomial Regression, Elastic Net, dan Generalized Additive Model (GAM) dalam memprediksi PE berdasarkan variabel AT, V, AP, dan RH pada dataset CCPP yang tersedia dalam beberapa sheet.
Adapun tujuan penelitian ini adalah: 1. Mendeskripsikan karakteristik data AT, V, AP, RH, dan PE pada setiap sheet. 2. Membandingkan kinerja Polynomial Regression, Elastic Net, dan GAM dalam memprediksi PE. 3. Memilih model terbaik berdasarkan AICc dan 10-fold Cross-Validation. 4. Menganalisis perbedaan kinerja model antar sheet menggunakan ANOVA dan Kruskal-Wallis. 5. Memberikan rekomendasi model prediksi PE yang paling sesuai berdasarkan hasil analisis.
1.4. Manfaat Penelitian
Secara praktis, hasil penelitian ini diharapkan dapat memberikan informasi mengenai model prediksi PE yang paling akurat dan stabil. Informasi tersebut dapat dimanfaatkan oleh operator CCPP atau pengelola data energi sebagai dasar pengambilan keputusan, optimasi operasional, dan evaluasi kinerja pembangkit.
Secara akademik, penelitian ini diharapkan memberikan kontribusi pada penerapan metode regresi dan pembelajaran mesin, khususnya perbandingan Polynomial Regression, Elastic Net, dan GAM dengan pemilihan model AICc dan 10-fold Cross-Validation. Penelitian ini juga dapat menjadi referensi untuk studi lanjutan pada dataset energi atau data dengan struktur multi-sheet.
1.5. Batasan Penelitian
Agar ruang lingkup penelitian lebih terarah, ditetapkan batasan sebagai berikut.
- Data yang digunakan adalah dataset regresi dengan variabel AT, V, AP, RH, dan PE.
- Dataset terdiri atas beberapa sheet (Sheet1–Sheet4) yang dianalisis secara terpisah.
- Metode yang dibandingkan adalah Polynomial Regression (derajat 1–3), Elastic Net, dan GAM.
- Pemilihan model menggunakan AICc dan 10-fold Cross-Validation.
- Evaluasi model uji menggunakan RMSE dan R².
- Uji beda kinerja antar sheet menggunakan ANOVA dan Kruskal-Wallis.
- Penelitian tidak membahas metode deep learning, neural network, atau pemodelan deret waktu.
Bab 2. Tinjauan Pustaka
2.1. Pendahuluan
Penelitian ini menggunakan pendekatan regresi untuk memprediksi PE berdasarkan AT, V, AP, dan RH. Pendekatan yang digunakan mencakup Polynomial Regression, Elastic Net, dan Generalized Additive Model (GAM). Bab ini membahas konsep dasar data mining, machine learning, regresi, regularisasi, GAM, pemilihan model, evaluasi kinerja, serta penelitian terdahulu pada dataset CCPP. Sesuai ketentuan, bab ini disajikan tanpa penulisan rumus.
2.2. Data Mining dan Machine Learning
Data mining adalah proses menggali pola, hubungan, atau pengetahuan baru dari data berukuran besar menggunakan pendekatan komputasional. Machine learning merupakan cabang kecerdasan buatan yang memungkinkan sistem belajar dari data. Dalam konteks prediksi, machine learning digunakan untuk membangun model yang mampu memetakan input ke output. Metode yang digunakan dalam penelitian ini termasuk dalam pendekatan supervised learning karena data memiliki variabel respons PE.
2.3. Regresi Linear dan Polynomial Regression
Regresi linear memodelkan hubungan antara variabel respons dan prediktor secara linear. Namun, hubungan antara kondisi lingkungan dan daya keluaran listrik pada CCPP sering tidak linear. Polynomial Regression memperluas regresi linear dengan menambahkan suku pangkat dan interaksi, sehingga mampu menangkap pola non-linear. Derajat polinomial menentukan fleksibilitas model; derajat terlalu rendah dapat menyebabkan underfitting, sedangkan derajat terlalu tinggi dapat menyebabkan overfitting.
2.4. Elastic Net
Elastic Net merupakan metode regularisasi yang menggabungkan penalti L1 (lasso) dan L2 (ridge). Penalti L1 mendorong seleksi variabel dengan menyusutkan koefisien tertentu menjadi nol, sedangkan penalti L2 menangani multikolinearitas dan menjaga stabilitas estimasi. Elastic Net sangat berguna ketika jumlah prediktor banyak atau terdapat korelasi antar prediktor. Parameter lambda mengatur besarnya penalti, sedangkan alpha mengatur proporsi antara L1 dan L2.
2.5. Generalized Additive Model (GAM)
Generalized Additive Model (GAM) merupakan perluasan model linear yang memodelkan hubungan antara respons dan prediktor melalui fungsi halus (smooth functions). GAM tidak mengasumsikan bentuk hubungan yang kaku, sehingga sangat fleksibel untuk menangkap non-linearitas. Dalam penelitian ini, GAM digunakan dengan fungsi spline dan parameter penghalus (smoothing parameter) yang dipilih melalui GCV atau grid search cross-validation.
2.6. Pemilihan Model: AICc dan Cross-Validation
Pemilihan model perlu mempertimbangkan keseimbangan antara kecocokan dan kompleksitas. AIC (Akaike Information Criterion) mengukur kecocokan model dengan memberikan penalti pada jumlah parameter. AICc merupakan koreksi AIC untuk ukuran sampel kecil. Semakin kecil nilai AICc, semakin baik model.
Cross-Validation (CV) mengukur kemampuan generalisasi model dengan membagi data menjadi beberapa fold. Pada 10-fold CV, data dibagi menjadi 10 bagian; model dilatih pada 9 bagian dan diuji pada 1 bagian secara bergantian. RMSE dari setiap fold dirata-ratakan untuk memperoleh ukuran kinerja yang lebih stabil.
2.7. Evaluasi Kinerja Model
Kinerja model prediksi dievaluasi menggunakan RMSE dan R². RMSE mengukur besarnya kesalahan prediksi dalam satuan variabel respons; semakin kecil RMSE, semakin baik. R² mengukur proporsi keragaman respons yang dapat dijelaskan oleh model; semakin mendekati 1, semakin baik.
2.8. Uji Beda Kinerja antar Sheet
Untuk mengetahui apakah kinerja model berbeda secara signifikan antar sheet, digunakan ANOVA jika asumsi normalitas dan homogenitas varians terpenuhi, atau Kruskal-Wallis sebagai alternatif non-parametrik. Uji ini diterapkan pada nilai CV-fold RMSE dari setiap sheet.
2.9. Penelitian Terdahulu pada Dataset CCPP
Dataset CCPP telah banyak digunakan dalam penelitian machine learning. Tüfekci (2014) membandingkan berbagai metode seperti regresi linear, neural network, dan support vector regression untuk memprediksi daya keluaran. Kaya et al. (2012) menggunakan metode local learning dan global learning. Penelitian-penelitian tersebut menunjukkan bahwa pendekatan non-linear dan regularisasi dapat meningkatkan akurasi prediksi. Namun, perbandingan khusus antara Polynomial Regression, Elastic Net, dan GAM dengan pemilihan model AICc dan 10-fold CV pada data multi-sheet masih terbatas.
2.10. Kerangka Pemikiran
Berdasarkan tinjauan pustaka, penelitian ini berangkat dari kebutuhan prediksi PE yang akurat. Tiga metode regresi dibandingkan pada setiap sheet. Model terbaik dipilih berdasarkan AICc dan 10-fold CV. Kinerja antar sheet dianalisis menggunakan ANOVA dan Kruskal-Wallis. Hasilnya digunakan untuk merekomendasikan model prediksi yang paling sesuai.
Bab 3. Metodologi Penelitian
3.1. Data Penelitian
Penelitian ini menggunakan data sekunder dari Combined Cycle Power Plant (CCPP) Data Set yang tersedia secara publik di UCI Machine Learning Repository. Dataset ini didonasikan pada 25 Maret 2014 dan pertama kali digunakan oleh Tüfekci (2014) untuk memprediksi daya keluaran listrik pada pembangkit siklus gabungan (Combined Cycle Power Plant). Dataset asli memuat 9.568 titik data yang dikumpulkan selama 6 tahun (2006–2011) pada kondisi pembangkit beroperasi dengan beban penuh (full load)
** Tabel 3.1. Variabel penelitian **
\[ \begin{array}{|l|l|l|} \hline \textbf{Variabel} & \textbf{Keterangan} & \textbf{Skala Data} \\ \hline AT & \text{Suhu lingkungan} & \text{Interval} \\ V & \text{Vakum buang} & \text{Interval} \\ AP & \text{Tekanan udara ambien} & \text{Interval} \\ RH & \text{Kelembapan relatif} & \text{Interval} \\ PE & \text{Daya keluaran listrik} & \text{Interval} \\ \hline \end{array} \] Variabel PE berperan sebagai variabel respons, sedangkan AT, V, AP, dan RH berperan sebagai prediktor. Data yang digunakan terdiri atas 5 dataset yang merupakan hasil partisi dari dataset CCPP asli ke dalam beberapa bagian untuk keperluan analisis multi-sheet dalam penelitian ini.
3.2 Metode Analisis
Penelitian ini menggunakan pendekatan kuantitatif untuk memodelkan Net Electrical Output (PE) berdasarkan variabel Ambient Temperature (AT), Exhaust Vacuum (V), Ambient Pressure (AP), dan Relative Humidity (RH). Analisis dilakukan pada setiap sheet menggunakan Polynomial Regression, Elastic Net, dan Generalized Additive Model (GAM). Pemilihan model dilakukan berdasarkan AICc dan 10-fold Cross-Validation, kemudian kinerja model dievaluasi menggunakan data uji.
3.2.1 Eksplorasi Data dan Pemeriksaan Awal
Eksplorasi data dilakukan untuk mengetahui karakteristik dan kualitas data. Analisis meliputi statistik deskriptif, pemeriksaan missing values, data duplikat, outlier, korelasi antarvariabel, dan multikolinearitas. Statistik deskriptif yang digunakan meliputi mean, standar deviasi, minimum, kuartil, median, maksimum, IQR, skewness, dan kurtosis.
Identifikasi outlier dilakukan menggunakan aturan Interquartile Range (IQR). Nilai IQR dihitung sebagai:
\[ IQR=Q_3-Q_1 \]
Batas bawah dan batas atas outlier ditentukan dengan:
\[ Batas\ bawah=Q_1-1.5(IQR) \]
\[ Batas\ atas=Q_3+1.5(IQR) \]
Skewness dan kurtosis digunakan untuk menggambarkan karakteristik distribusi data. Skewness dihitung dengan:
\[ Skewness= \frac{\frac{1}{n}\sum_{i=1}^{n}(x_i-\bar{x})^3} {\left(\frac{1}{n}\sum_{i=1}^{n}(x_i-\bar{x})^2\right)^{3/2}} \]
Kurtosis dihitung dengan:
\[ Kurtosis= \frac{\frac{1}{n}\sum_{i=1}^{n}(x_i-\bar{x})^4} {\left(\frac{1}{n}\sum_{i=1}^{n}(x_i-\bar{x})^2\right)^2}-3 \]
Multikolinearitas antarvariabel prediktor diperiksa menggunakan Variance Inflation Factor (VIF) dengan:
\[ VIF_j=\frac{1}{1-R_j^2} \]
dengan \(R_j^2\) merupakan koefisien determinasi dari regresi prediktor ke-\(j\) terhadap prediktor lainnya.
3.2.2 Pembagian Data dan 10-Fold Cross-Validation
Data pada setiap sheet dibagi menjadi 80% data latih dan 20% data uji secara acak menggunakan set.seed(42). Jika \(n\) merupakan jumlah observasi, jumlah data uji ditentukan dengan:
\[ n_{test}=\lceil0.20n\rceil \]
sedangkan jumlah data latih adalah:
\[ n_{train}=n-n_{test} \]
Data latih digunakan untuk membangun model dan melakukan 10-fold Cross-Validation, sedangkan data uji digunakan untuk mengevaluasi kinerja akhir model. Pada setiap iterasi Cross-Validation, sembilan fold digunakan sebagai data pelatihan dan satu fold sebagai data validasi.
3.2.3 Polynomial Regression
Polynomial Regression digunakan untuk menangkap hubungan nonlinier antara variabel prediktor dan respons. Model polynomial dengan \(p\) prediktor dan derajat \(d\) dinyatakan sebagai:
\[ y_i=\beta_0+ \sum_{j=1}^{p}\sum_{r=1}^{d} \beta_{jr}x_{ij}^{r}+\varepsilon_i \]
Dalam penelitian ini digunakan derajat polynomial 1, 2, dan 3. Estimasi parameter dilakukan menggunakan Ordinary Least Squares (OLS) dengan:
\[ \hat{\beta}=(X^TX)^{-1}X^Ty \]
Residual Sum of Squares (RSS) digunakan sebagai ukuran kesalahan model:
\[ RSS=\sum_{i=1}^{n}(y_i-\hat{y}_i)^2 \]
3.2.4 Elastic Net
Elastic Net digunakan untuk mengatasi kompleksitas model dan potensi multikolinearitas dengan menggabungkan penalti L1 dan L2. Fungsi objektifnya adalah:
\[ \min_{\beta_0,\beta} \left\{ \frac{1}{2n}\sum_{i=1}^{n} (y_i-\beta_0-x_i^T\beta)^2 + \lambda \left[ \frac{1-\alpha}{2}\|\beta\|_2^2+ \alpha\|\beta\|_1 \right] \right\} \]
dengan \(\lambda\) sebagai parameter penalti dan \(\alpha\) sebagai proporsi penalti L1. Nilai \(\alpha=1\) menghasilkan Lasso, sedangkan \(\alpha=0\) menghasilkan Ridge.
Sebelum pemodelan, variabel prediktor distandarisasi berdasarkan data latih dengan:
\[ x_{scaled}=\frac{x-\bar{x}}{s} \]
Parameter \(\lambda\) dan \(\alpha\) ditentukan melalui kombinasi nilai yang diuji pada proses pemodelan.
3.2.5 Generalized Additive Model (GAM)
Generalized Additive Model (GAM) digunakan untuk menangkap hubungan nonlinier melalui fungsi smooth. Bentuk umum model GAM adalah:
\[ y_i=\beta_0+\sum_{j=1}^{p}f_j(x_{ij})+\varepsilon_i \]
Dalam penelitian ini, fungsi smooth digunakan pada masing-masing prediktor dengan basis spline:
\[ y\sim s(AT,k=20)+s(V,k=20)+s(AP,k=20)+s(RH,k=20) \]
Fungsi smooth dibentuk dari kombinasi fungsi basis:
\[ f_j(x)=\sum_{m=1}^{k}\beta_{jm}b_{jm}(x) \]
Kehalusan fungsi dikendalikan menggunakan parameter penghalus \(\lambda\). Estimasi penghalusan juga mempertimbangkan Generalized Cross-Validation (GCV):
\[ GCV= \frac{RSS/n}{(1-edf/n)^2} \]
dengan \(edf\) merupakan effective degrees of freedom.
3.2.6 Pemilihan Model dengan AICc
Pemilihan model berdasarkan Akaike Information Criterion corrected (AICc) dilakukan dengan mempertimbangkan kesesuaian dan kompleksitas model. Nilai AIC dihitung dengan:
\[ AIC=n\log\left(\frac{RSS}{n}\right)+2k \]
Kemudian AICc dihitung sebagai:
\[ AICc=AIC+\frac{2k(k+1)}{n-k-1} \]
dengan \(n\) merupakan jumlah observasi dan \(k\) merupakan jumlah parameter model. Model dengan nilai AICc yang lebih kecil dipilih berdasarkan kriteria AICc.
3.2.7 Evaluasi 10-Fold Cross-Validation
Kinerja model dalam proses Cross-Validation dievaluasi menggunakan Root Mean Squared Error (RMSE). RMSE pada fold ke-\(k\) dihitung dengan:
\[ RMSE_k= \sqrt{ \frac{1}{n_k} \sum_{i\in fold_k} (y_i-\hat{y}_i)^2 } \]
Nilai CV_RMSE diperoleh dari akar rata-rata MSE seluruh fold:
\[ CV\_RMSE= \sqrt{ \frac{1}{K} \sum_{k=1}^{K}MSE_k } \]
dengan \(K=10\). Model dengan nilai CV_RMSE yang lebih kecil menunjukkan kesalahan prediksi yang lebih rendah.
3.2.8 Evaluasi Model pada Data Uji
Model terpilih selanjutnya dievaluasi menggunakan data uji berdasarkan RMSE dan R². RMSE dihitung dengan:
\[ RMSE= \sqrt{ \frac{1}{n_{test}} \sum_{i=1}^{n_{test}} (y_i-\hat{y}_i)^2 } \]
Sedangkan koefisien determinasi dihitung dengan:
\[ R^2= 1- \frac{ \sum_{i=1}^{n_{test}}(y_i-\hat{y}_i)^2 }{ \sum_{i=1}^{n_{test}}(y_i-\bar{y})^2 } \]
RMSE menunjukkan besarnya kesalahan prediksi, sedangkan R² menunjukkan proporsi variasi PE yang dapat dijelaskan oleh model.
3.2.9 Korelasi Peringkat Spearman
Korelasi peringkat Spearman digunakan untuk mengetahui kesesuaian peringkat model berdasarkan AICc dan CV_RMSE. Koefisien Spearman dihitung dengan:
\[ \rho= 1-\frac{6\sum_{i=1}^{n}d_i^2} {n(n^2-1)} \]
dengan \(d_i\) merupakan selisih peringkat antara kedua kriteria dan \(n\) merupakan jumlah model yang dibandingkan.
3.2.10 ANOVA dan Kruskal-Wallis
Perbedaan kinerja model antar-sheet dianalisis menggunakan ANOVA satu arah dan Kruskal-Wallis berdasarkan nilai CV-fold RMSE.
Statistik ANOVA dihitung dengan:
\[ F= \frac{SSB/(k-1)} {SSW/(N-k)} \]
dengan \(SSB\) merupakan jumlah kuadrat antar kelompok dan \(SSW\) merupakan jumlah kuadrat dalam kelompok.
Apabila diperlukan pendekatan nonparametrik, digunakan uji Kruskal-Wallis dengan statistik:
\[ H= \frac{12}{N(N+1)} \sum_{i=1}^{k}\frac{R_i^2}{n_i} -3(N+1) \]
Kedua pengujian digunakan untuk mengetahui apakah terdapat perbedaan kinerja model antar-sheet.
3.2.11 Uji Normalitas dan Heteroskedastisitas
Pemeriksaan asumsi dilakukan menggunakan uji Shapiro-Wilk untuk normalitas residual dan uji Breusch-Pagan untuk heteroskedastisitas.
Statistik Shapiro-Wilk dinyatakan sebagai:
\[ W= \frac{ \left(\sum_{i=1}^{n}a_ix_{(i)}\right)^2 }{ \sum_{i=1}^{n}(x_i-\bar{x})^2 } \]
Sementara itu, statistik Breusch-Pagan dihitung dengan:
\[ LM=nR^2 \]
Hasil kedua pengujian digunakan sebagai informasi pendukung dalam interpretasi model.
3.2.12 Tahapan Penelitian
Tahapan penelitian dilakukan secara sistematis melalui proses
- membaca dan memvalidasi data,
- memeriksa kesamaan data antar-sheet,
- melakukan eksplorasi data,
- membagi data menjadi data latih dan data uji,
- melakukan 10-fold Cross-Validation,
- membangun Polynomial Regression, Elastic Net, dan GAM,
- melakukan pemilihan model berdasarkan AICc dan CV_RMSE,
- mengevaluasi model menggunakan RMSE dan R² pada data uji,
- menganalisis kesesuaian peringkat dan perbedaan kinerja antar-sheet menggunakan korelasi Spearman, ANOVA, dan Kruskal-Wallis.
Bab 4. Hasil Penelitian
4.1. Hasil Eksplorasi dan Pemeriksaan Data
Tahap eksplorasi dan pemeriksaan data dilakukan untuk mengetahui struktur, karakteristik, serta kondisi awal data sebelum proses pemodelan. Analisis dilakukan secara berurutan meliputi pembacaan dan validasi data, pemeriksaan kesamaan antar-sheet, statistik deskriptif, pemeriksaan missing value dan duplikasi, identifikasi outlier, analisis hubungan antarvariabel, pemeriksaan multikolinearitas, serta pemeriksaan residual pada model regresi awal.
4.1.1. Struktur dan Karakteristik Data
Hasil pemeriksaan menunjukkan bahwa seluruh data memiliki struktur variabel yang sama dan jumlah observasi yang konsisten. Pemeriksaan kesamaan data setelah pengurutan juga menunjukkan hasil yang identik. Statistik deskriptif menunjukkan bahwa PE memiliki rata-rata sebesar 454,365 dengan standar deviasi 17,067. Variabel AT, V, AP, dan RH memiliki tingkat penyebaran yang berbeda sesuai karakteristik masing-masing. Nilai skewness seluruh variabel relatif kecil, sehingga distribusi data tidak menunjukkan kemencengan yang sangat kuat.
Syntax yang digunakan:
fname <- "DATASET REGRESI.xlsx"
sheets <- readxl::excel_sheets(fname)
data_sheets <- lapply(sheets, function(sh) {
d <- readxl::read_excel(fname, sheet = sh)
colnames(d) <- trimws(as.character(colnames(d)))
as.data.frame(d)
})
names(data_sheets) <- sheets
sheet_names <- names(data_sheets)
target <- "PE"
features <- c("AT", "V", "AP", "RH")
validate_ccpp_data <- function(data_sheets, target, features) {
if (length(data_sheets) == 0) {
stop("Tidak ada sheet yang terbaca.")
}
need <- c(features, target)
problems <- c()
for (sh in names(data_sheets)) {
d <- data_sheets[[sh]]
missing <- setdiff(need, colnames(d))
if (length(missing) > 0) {
problems <- c(
problems,
paste0(sh, ": ", paste(missing, collapse = ", "))
)
}
}
if (length(problems) > 0) {
stop(
paste(
"Variabel tidak ditemukan ->",
paste(problems, collapse = " | ")
)
)
}
}
validate_ccpp_data(data_sheets, target, features)
df <- data_sheets[[1]]verify_identical_sheets <- function(data_sheets, features, target) {
cols <- c(features, target)
sheet_names <- names(data_sheets)
ref <- data_sheets[[sheet_names[1]]][, cols]
ref <- ref[do.call(order, ref), ]
rownames(ref) <- NULL
rows <- list()
for (sh in sheet_names) {
d_sorted <- data_sheets[[sh]][, cols]
d_sorted <- d_sorted[do.call(order, d_sorted), ]
rownames(d_sorted) <- NULL
identik <- identical(ref, d_sorted)
rows[[sh]] <- data.frame(
Sheet = sh,
n_baris = nrow(data_sheets[[sh]]),
identik_dgn_sheet1_setelah_sort = identik
)
}
return(do.call(rbind, rows))
}
verify_df <- verify_identical_sheets(data_sheets, features, target)
verify_df## Sheet n_baris identik_dgn_sheet1_setelah_sort
## Sheet1 Sheet1 9568 TRUE
## Sheet2 Sheet2 9568 TRUE
## Sheet3 Sheet3 9568 TRUE
## Sheet4 Sheet4 9568 TRUE
## Sheet5 Sheet5 9568 TRUE# Ukuran pusat & sebaran lengkap
calc_mode <- function(x) {
ux <- unique(x)
ux[which.max(tabulate(match(x, ux)))]
}
calc_skewness <- function(x) {
n <- length(x)
m3 <- sum((x - mean(x))^3) / n
s3 <- (sum((x - mean(x))^2) / n)^(3/2)
return(m3 / s3)
}
calc_kurtosis <- function(x) {
n <- length(x)
m4 <- sum((x - mean(x))^4) / n
s4 <- (sum((x - mean(x))^2) / n)^2
return((m4 / s4) - 3)
}
desc_stats <- data.frame(
Mean = colMeans(df),
Std = apply(df, 2, sd),
Min = apply(df, 2, min),
Q25 = apply(df, 2, quantile, probs = 0.25),
Median = apply(df, 2, median),
Mode = apply(df, 2, calc_mode),
Q75 = apply(df, 2, quantile, probs = 0.75),
Max = apply(df, 2, max),
IQR = apply(df, 2, IQR),
Skewness = apply(df, 2, calc_skewness),
Kurtosis = apply(df, 2, calc_kurtosis)
)
desc_stats## Mean Std Min Q25 Median Mode Q75 Max
## AT 19.65123 7.452473 1.81 13.5100 20.345 25.21 25.72 37.11
## V 54.30580 12.707893 25.36 41.7400 52.080 41.17 66.54 81.56
## AP 1013.25908 5.938784 992.89 1009.1000 1012.940 1013.88 1017.26 1033.30
## RH 73.30898 14.600269 25.56 63.3275 74.975 100.09 84.83 100.16
## PE 454.36501 17.066995 420.26 439.7500 451.550 468.80 468.43 495.76
## IQR Skewness Kurtosis
## AT 12.2100 -0.1363717 -1.03763410
## V 24.8000 0.1984899 -1.44420904
## AP 8.1600 0.2654031 0.09356093
## RH 21.5025 -0.4317710 -0.44492113
## PE 28.6800 0.3064614 -1.048600144.1.2. Pemeriksaan Kualitas dan Distribusi Data
Hasil pemeriksaan menunjukkan tidak terdapat missing value. Ditemukan 41 observasi duplikat. Berdasarkan aturan IQR, outlier teridentifikasi pada AP sebanyak 88 observasi dan RH sebanyak 12 observasi, sedangkan AT, V, dan PE tidak menunjukkan outlier berdasarkan kriteria tersebut.
Histogram dan boxplot digunakan untuk melihat bentuk distribusi dan penyebaran data secara visual. Hasil visualisasi menunjukkan bahwa setiap variabel memiliki pola distribusi dan tingkat penyebaran yang berbeda, dengan beberapa nilai ekstrem terutama pada AP dan RH.
Syntax yang digunakan:
colSums(is.na(df))## AT V AP RH PE
## 0 0 0 0 0cat("\nJumlah baris duplikat:", sum(duplicated(df)), "\n")##
## Jumlah baris duplikat: 41p_list <- lapply(colnames(df), function(col) {
ggplot(df, aes_string(x = col)) +
geom_histogram(aes(y = ..density..), fill = "steelblue", bins = 30, alpha = 0.7) +
geom_density(color = "darkblue", size = 0.8) +
labs(title = paste("Distribusi", col)) +
theme_minimal()
})
grid.arrange(grobs = p_list, ncol = 3)
p_box <- lapply(colnames(df), function(col) {
ggplot(df, aes_string(y = col)) +
geom_boxplot(fill = "lightcoral") +
labs(title = col) +
theme_minimal()
})
grid.arrange(grobs = p_box, ncol = 5)
iqr_outlier_count <- function(s) {
q1 <- quantile(s, 0.25)
q3 <- quantile(s, 0.75)
iqr <- q3 - q1
lb <- q1 - 1.5 * iqr
ub <- q3 + 1.5 * iqr
sum(s < lb | s > ub)
}
print(sapply(df, iqr_outlier_count))## AT V AP RH PE
## 0 0 88 12 04.1.3. Hubungan Antarvariabel dan Pemeriksaan Asumsi Awal
Hasil analisis korelasi menunjukkan adanya hubungan yang kuat antara beberapa variabel prediktor dengan PE. AT dan V memiliki hubungan negatif yang kuat dengan PE, sedangkan AP dan RH menunjukkan hubungan positif dengan PE. Hubungan yang cukup kuat juga ditemukan antara AT dan V.
Pemeriksaan multikolinearitas menghasilkan nilai VIF sebesar 5,978 untuk AT, 3,943 untuk V, 1,453 untuk AP, dan 1,705 untuk RH. Nilai tersebut menunjukkan bahwa AT memiliki tingkat keterkaitan dengan prediktor lainnya paling tinggi dibandingkan variabel prediktor lainnya.
Pemeriksaan residual pada model regresi awal menghasilkan statistik Shapiro-Wilk sebesar 0,9795 dengan p-value sebesar 2,728 × 10⁻²⁶ dan uji Breusch-Pagan sebesar 33,9302 dengan p-value sebesar 7,702 × 10⁻⁷. Hasil tersebut menunjukkan adanya penyimpangan dari normalitas residual dan indikasi heteroskedastisitas. Temuan ini menjadi pertimbangan untuk menggunakan metode pemodelan yang mampu menangkap hubungan nonlinier serta mengakomodasi kompleksitas hubungan antarvariabel.
Syntax yang digunakan:
cor_matrix <- cor(df)
corrplot(cor_matrix, method = "color", addCoef.col = "black",
tl.col = "black", number.digits = 2, cl.lim = c(-1, 1),
title = "Matriks Korelasi", mar = c(0,0,2,0))
p_reg <- lapply(features, function(col) {
ggplot(df, aes_string(x = col, y = target)) +
geom_point(alpha = 0.15, size = 1) +
geom_smooth(method = "lm", color = "red", se = FALSE) +
labs(title = paste(col, "vs", target)) +
theme_minimal()
})
grid.arrange(grobs = p_reg, ncol = 4)
X_check <- lm(as.formula(paste(target, "~", paste(features, collapse = " + "))), data = df)
vif_vals <- car::vif(X_check)
vif_df <- data.frame(
variabel = names(vif_vals),
VIF = as.numeric(vif_vals)
)
vif_df## variabel VIF
## 1 AT 5.977602
## 2 V 3.943003
## 3 AP 1.452639
## 4 RH 1.705290model0 <- lm(as.formula(paste(target, "~", paste(features, collapse = " + "))), data = df)
resid <- residuals(model0)
set.seed(RANDOM_STATE)
sample_resid <- sample(resid, min(5000, length(resid)))
shapiro_res <- shapiro.test(sample_resid)
cat(sprintf("Shapiro-Wilk: statistic=%.4f, p-value=%.4g\n", shapiro_res$statistic, shapiro_res$p.value))## Shapiro-Wilk: statistic=0.9795, p-value=2.728e-26p1 <- ggplot(data.frame(resid = resid), aes(x = resid)) +
geom_histogram(aes(y = ..density..), bins = 30, fill = "steelblue", alpha = 0.7) +
geom_density(color = "darkblue") +
labs(title = "Distribusi Residual") + theme_minimal()
p2 <- ggplot(data.frame(resid = resid), aes(sample = resid)) +
stat_qq() + stat_qq_line(color = "red") +
labs(title = "QQ-Plot Residual") + theme_minimal()
grid.arrange(p1, p2, ncol = 2)
bp_res <- lmtest::bptest(model0)
cat(sprintf("Breusch-Pagan: statistic=%.4f, p-value=%.4g\n", bp_res$statistic, bp_res$p.value))## Breusch-Pagan: statistic=33.9302, p-value=7.702e-07cat("(p-value kecil -> indikasi heteroskedastisitas, perlu dicatat sebagai catatan interpretasi, bukan alasan transformasi otomatis)\n")## (p-value kecil -> indikasi heteroskedastisitas, perlu dicatat sebagai catatan interpretasi, bukan alasan transformasi otomatis)compute_aic <- function(n, rss, k) {
n * log(rss / n) + 2 * k
}
compute_aicc <- function(n, rss, k) {
aic <- compute_aic(n, rss, k)
if (n - k - 1 > 0) {
aic <- aic + (2 * k * (k + 1)) / (n - k - 1)
}
return(aic)
}
rmse <- function(y_true, y_pred) {
sqrt(mean((y_true - y_pred)^2))
}
r2_score <- function(y_true, y_pred) {
1 - sum((y_true - y_pred)^2) / sum((y_true - mean(y_true))^2)
}
make_train_test_split <- function(n, test_size = 0.20, seed = 42) {
set.seed(seed)
shuffled <- sample(n)
n_test <- ceiling(n * test_size)
test_idx <- shuffled[1:n_test]
train_idx <- shuffled[(n_test + 1):n]
return(list(train = train_idx, test = test_idx))
}
make_folds <- function(n, K = 10, seed = 42) {
set.seed(seed)
shuffled <- sample(n)
fold_id <- rep(1:K, length.out = n)
folds <- lapply(1:K, function(k) sort(shuffled[fold_id == k]))
return(folds)
}
fit_scaler <- function(X) {
center <- colMeans(X)
scale <- apply(X, 2, sd)
scale[scale == 0 | is.na(scale)] <- 1.0
return(list(center = center, scale = scale))
}
apply_scaler <- function(X, center, scale) {
scale(X, center = center, scale = scale)
}4.2. Hasil Pemodelan
Setelah seluruh model dijalankan, hasil pemodelan diperoleh untuk setiap sheet dengan menggunakan tiga metode, yaitu Polynomial Regression, Elastic Net, dan Generalized Additive Model (GAM). Proses pemodelan dilakukan menggunakan pembagian data dan skema 10-fold Cross-Validation yang sama sehingga hasil antar-metode dapat dibandingkan.
Hasil pemodelan menunjukkan bahwa ketiga metode menghasilkan nilai kesalahan prediksi yang berbeda. Nilai CV_RMSE Polynomial berada pada kisaran sekitar 4,23–4,28, sedangkan Elastic Net berada pada kisaran sekitar 4,57–4,58. GAM menghasilkan CV_RMSE sekitar 4,12–4,16. Berdasarkan nilai tersebut, GAM menghasilkan nilai CV_RMSE yang lebih rendah dibandingkan Polynomial Regression dan Elastic Net pada seluruh sheet.
Untuk melihat konfigurasi dan nilai kriteria masing-masing metode secara lebih lengkap, hasil pemodelan selanjutnya disusun berdasarkan AICc dan CV_RMSE. Hasil tersebut digunakan pada tahap perbandingan dan pemilihan model pada subbab berikutnya.
Syntax yang digunakan:
run_polynomial <- function(Xtr, ytr, Xte, yte, folds) {
n <- length(ytr)
rows <- list()
fold_scores <- list()
# ============================================================
# Polynomial Regression degree 1-3
# ============================================================
for (degree in 1:3) {
# Membuat polynomial features secara eksplisit
Xtr_df <- as.data.frame(Xtr)
Xte_df <- as.data.frame(Xte)
# Formula polynomial
poly_terms <- unlist(
lapply(colnames(Xtr_df), function(v) {
paste0(
"poly(", v, ", degree = ", degree,
", raw = TRUE)"
)
})
)
formula_poly <- as.formula(
paste("~", paste(poly_terms, collapse = " + "))
)
# Design matrix training dan testing
Xp_tr <- model.matrix(formula_poly, data = Xtr_df)
Xp_te <- model.matrix(formula_poly, data = Xte_df)
# ------------------------------------------------------------
# Model pada seluruh training data
# ------------------------------------------------------------
lr <- lm(ytr ~ Xp_tr - 1)
rss <- sum(residuals(lr)^2)
# Jumlah parameter
k <- ncol(Xp_tr) + 1
aicc <- compute_aicc(
n = n,
rss = rss,
k = k
)
# ------------------------------------------------------------
# 10-fold Cross Validation
# ------------------------------------------------------------
fold_rmse <- numeric(length(folds))
for (i in seq_along(folds)) {
val_idx <- folds[[i]]
tr_idx <- setdiff(1:n, val_idx)
fold_fit <- lm(
ytr[tr_idx] ~ Xp_tr[tr_idx, , drop = FALSE] - 1
)
pred_val <- Xp_tr[val_idx, , drop = FALSE] %*%
coef(fold_fit)
pred_val <- as.numeric(pred_val)
# Penanganan koefisien NA
pred_val[is.na(pred_val)] <- mean(ytr[tr_idx])
fold_rmse[i] <- rmse(
ytr[val_idx],
pred_val
)
}
# Simpan hasil CV setiap degree
fold_scores[[as.character(degree)]] <- fold_rmse
rows[[degree]] <- data.frame(
degree = degree,
k = k,
AICc = aicc,
CV_RMSE = sqrt(mean(fold_rmse^2))
)
}
# ============================================================
# Tabel hasil pemilihan degree
# ============================================================
res <- do.call(rbind, rows)
degree_aic <- res$degree[
which.min(res$AICc)
]
degree_cv <- res$degree[
which.min(res$CV_RMSE)
]
# ============================================================
# Evaluasi model terpilih
# ============================================================
evaluate_route <- function(route, degree) {
Xtr_df <- as.data.frame(Xtr)
Xte_df <- as.data.frame(Xte)
poly_terms <- unlist(
lapply(colnames(Xtr_df), function(v) {
paste0(
"poly(", v, ", degree = ", degree,
", raw = TRUE)"
)
})
)
formula_poly <- as.formula(
paste("~", paste(poly_terms, collapse = " + "))
)
Xp_tr <- model.matrix(
formula_poly,
data = Xtr_df
)
Xp_te <- model.matrix(
formula_poly,
data = Xte_df
)
fit <- lm(
ytr ~ Xp_tr - 1
)
pred_test <- Xp_te %*%
coef(fit)
pred_test <- as.numeric(pred_test)
pred_test[is.na(pred_test)] <- mean(ytr)
row <- res[
res$degree == degree,
,
drop = FALSE
]
list(
config = paste0(
"degree=", degree
),
AICc = if (
route == "AIC"
) row$AICc else NA,
CV_RMSE = if (
route == "CV"
) row$CV_RMSE else NA,
test_RMSE = rmse(
yte,
pred_test
),
test_R2 = r2_score(
yte,
pred_test
)
)
}
# ============================================================
# Output
# ============================================================
list(
table = res,
final = list(
AIC = evaluate_route(
"AIC",
degree_aic
),
CV = evaluate_route(
"CV",
degree_cv
)
),
cv_fold_rmse = fold_scores
)
}run_elasticnet <- function(Xtr, ytr, Xte, yte, folds,
lambda_grid = NULL, l1_ratio_grid = c(0.1, 0.5, 0.9, 1.0), seed = 42) {
if (is.null(lambda_grid)) {
lambda_grid <- sort(10^seq(-4, 1, length.out = 20), decreasing = TRUE)
}
scaler <- fit_scaler(Xtr)
Xtr_s <- apply_scaler(Xtr, scaler$center, scaler$scale)
Xte_s <- apply_scaler(Xte, scaler$center, scaler$scale)
n <- nrow(Xtr_s)
# 1. Grid Search AICc
aic_rows <- list()
for (l1r in l1_ratio_grid) {
for (lam in lambda_grid) {
m <- glmnet(Xtr_s, ytr, alpha = l1r, lambda = lam, standardise = FALSE)
pred <- predict(m, newx = Xtr_s)
rss <- sum((ytr - pred)^2)
k <- sum(as.vector(coef(m)) != 0)
aic_rows[[length(aic_rows) + 1]] <- data.frame(
l1_ratio = l1r, lambda = lam, k = k, AICc = compute_aicc(n, rss, k)
)
}
}
en_aic <- do.call(rbind, aic_rows)
best_aicc <- en_aic[which.min(en_aic$AICc), ]
# 2. Grid Search Cross-Validation
cv_rows <- list()
for (l1r in l1_ratio_grid) {
for (lam in lambda_grid) {
fold_mse <- numeric(length(folds))
for (i in seq_along(folds)) {
val_idx <- folds[[i]]
tr_idx <- setdiff(1:n, val_idx)
m <- glmnet(Xtr_s[tr_idx, ], ytr[tr_idx], alpha = l1r, lambda = lam, standardise = FALSE)
pred <- predict(m, newx = Xtr_s[val_idx, ])
fold_mse[i] <- mean((ytr[val_idx] - pred)^2)
}
cv_rows[[length(cv_rows) + 1]] <- data.frame(
l1_ratio = l1r, lambda = lam, CV_RMSE = sqrt(mean(fold_mse))
)
}
}
cv_grid <- do.call(rbind, cv_rows)
best_cv <- cv_grid[which.min(cv_grid$CV_RMSE), ]
# Hitung RMSE per fold untuk model CV terbaik
cv_fold_rmse <- numeric(length(folds))
for (i in seq_along(folds)) {
val_idx <- folds[[i]]
tr_idx <- setdiff(1:n, val_idx)
m <- glmnet(Xtr_s[tr_idx, ], ytr[tr_idx], alpha = best_cv$l1_ratio, lambda = best_cv$lambda, standardise = FALSE)
pred <- predict(m, newx = Xtr_s[val_idx, ])
cv_fold_rmse[i] <- rmse(ytr[val_idx], pred)
}
# 3. Evaluasi Test Set
fit_aic <- glmnet(Xtr_s, ytr, alpha = best_aicc$l1_ratio, lambda = best_aicc$lambda, standardise = FALSE)
pred_aic <- predict(fit_aic, newx = Xte_s)
fit_cv <- glmnet(Xtr_s, ytr, alpha = best_cv$l1_ratio, lambda = best_cv$lambda, standardise = FALSE)
pred_cv <- predict(fit_cv, newx = Xte_s)
list(
table = en_aic,
cv_grid = cv_grid,
final = list(
AIC = list(
config = sprintf("lambda=%.4f, l1=%.1f", best_aicc$lambda, best_aicc$l1_ratio),
AICc = best_aicc$AICc, CV_RMSE = NA,
test_RMSE = rmse(yte, pred_aic), test_R2 = r2_score(yte, pred_aic)
),
CV = list(
config = sprintf("lambda=%.4f, l1=%.1f", best_cv$lambda, best_cv$l1_ratio),
AICc = NA, CV_RMSE = mean(cv_fold_rmse),
test_RMSE = rmse(yte, pred_cv), test_R2 = r2_score(yte, pred_cv)
)
),
cv_fold_rmse = cv_fold_rmse,
scaler = scaler
)
}run_gam <- function(Xtr, ytr, Xte, yte, folds,
lambda_grid = NULL, lambda_grid_gcv = NULL, k_spline = 20) {
if (is.null(lambda_grid)) lambda_grid <- 10^seq(-3, 3, length.out = 10)
df_tr <- data.frame(y = ytr, Xtr)
df_te <- data.frame(y = yte, Xte)
gam_formula <- as.formula(paste("y ~", paste(sprintf("s(%s, k=%d)", colnames(Xtr), k_spline), collapse = " + ")))
gam_gcv <- gam(gam_formula, data = df_tr, method = "GCV.Cp")
rss <- sum(residuals(gam_gcv)^2)
edof <- sum(influence(gam_gcv))
aicc_gcv <- compute_aicc(nrow(df_tr), rss, edof)
lambda_gcv <- mean(gam_gcv$sp)
cv_rows <- list()
fold_scores <- list()
for (lam in lambda_grid) {
fold_rmse <- numeric(length(folds))
sp_vec <- rep(lam, length(colnames(Xtr)))
for (i in seq_along(folds)) {
val_idx <- folds[[i]]
tr_idx <- setdiff(1:nrow(df_tr), val_idx)
fit <- gam(gam_formula, data = df_tr[tr_idx, ], sp = sp_vec)
pred <- predict(fit, newdata = df_tr[val_idx, ])
fold_rmse[i] <- rmse(df_tr$y[val_idx], pred)
}
fold_scores[[as.character(lam)]] <- fold_rmse
cv_rows[[length(cv_rows) + 1]] <- data.frame(lambda = lam, CV_RMSE = mean(fold_rmse))
}
cv_df <- do.call(rbind, cv_rows)
best_cv <- cv_df[which.min(cv_df$CV_RMSE), ]
best_cv_fold_rmse <- fold_scores[[as.character(best_cv$lambda)]]
pred_test_gcv <- predict(gam_gcv, newdata = df_te)
gam_cv <- gam(gam_formula, data = df_tr, sp = rep(best_cv$lambda, length(colnames(Xtr))))
pred_test_cv <- predict(gam_cv, newdata = df_te)
list(
table = cv_df,
final = list(
AIC = list(config = sprintf("lambda(GCV)~=%.4g", lambda_gcv), AICc = aicc_gcv, CV_RMSE = NA,
test_RMSE = rmse(yte, pred_test_gcv), test_R2 = r2_score(yte, pred_test_gcv)),
CV = list(config = sprintf("lambda=%.4g", best_cv$lambda), AICc = NA, CV_RMSE = best_cv$CV_RMSE,
test_RMSE = rmse(yte, pred_test_cv), test_R2 = r2_score(yte, pred_test_cv))
),
cv_fold_rmse = best_cv_fold_rmse,
gcv_model = gam_gcv,
cv_model = gam_cv
)
}run_full_analysis <- function(data_sheets, target = "PE", features = c("AT", "V", "AP", "RH"), K = 10, seed = 42) {
sheet_names <- names(data_sheets)
n <- nrow(data_sheets[[sheet_names[1]]])
split_res <- make_train_test_split(n, test_size = 0.20, seed = seed)
train_idx <- split_res$train
test_idx <- split_res$test
folds <- make_folds(length(train_idx), K = K, seed = seed)
splits <- list()
for (sh in sheet_names) {
d <- data_sheets[[sh]]
splits[[sh]] <- list(
X_train = d[train_idx, features],
X_test = d[test_idx, features],
y_train = d[train_idx, target],
y_test = d[test_idx, target]
)
}
poly_all <- list()
en_all <- list()
gam_all <- list()
for (sh in sheet_names) {
sp <- splits[[sh]]
cat(sprintf("[%s] menjalankan Polynomial...\n", sh))
poly_all[[sh]] <- run_polynomial(sp$X_train, sp$y_train, sp$X_test, sp$y_test, folds)
cat(sprintf("[%s] menjalankan Elastic Net...\n", sh))
en_all[[sh]] <- run_elasticnet(sp$X_train, sp$y_train, sp$X_test, sp$y_test, folds, seed = seed)
cat(sprintf("[%s] menjalankan GAM...\n", sh))
gam_all[[sh]] <- run_gam(sp$X_train, sp$y_train, sp$X_test, sp$y_test, folds)
}
list(sheet_names = sheet_names, splits = splits, folds = folds,
poly_all = poly_all, en_all = en_all, gam_all = gam_all)
}
result <- run_full_analysis(data_sheets, target = target, features = features, K = 10, seed = RANDOM_STATE)## [Sheet1] menjalankan Polynomial...
## [Sheet1] menjalankan Elastic Net...
## [Sheet1] menjalankan GAM...
## [Sheet2] menjalankan Polynomial...
## [Sheet2] menjalankan Elastic Net...
## [Sheet2] menjalankan GAM...
## [Sheet3] menjalankan Polynomial...
## [Sheet3] menjalankan Elastic Net...
## [Sheet3] menjalankan GAM...
## [Sheet4] menjalankan Polynomial...
## [Sheet4] menjalankan Elastic Net...
## [Sheet4] menjalankan GAM...
## [Sheet5] menjalankan Polynomial...
## [Sheet5] menjalankan Elastic Net...
## [Sheet5] menjalankan GAM...cat("\nSelesai menjalankan Polynomial/Elastic Net/GAM untuk seluruh sheet:", paste(result$sheet_names, collapse = ", "), "\n")##
## Selesai menjalankan Polynomial/Elastic Net/GAM untuk seluruh sheet: Sheet1, Sheet2, Sheet3, Sheet4, Sheet5get_config_degree <- function(config) {
as.integer(sub(".*=", "", config))
}
build_best_model_df <- function(result) {
rows <- list()
for (sh in result$sheet_names) {
aicc_vals <- c(
Polynomial = result$poly_all[[sh]]$final$AIC$AICc,
"Elastic Net" = result$en_all[[sh]]$final$AIC$AICc,
GAM = result$gam_all[[sh]]$final$AIC$AICc
)
cv_vals <- c(
Polynomial = result$poly_all[[sh]]$final$CV$CV_RMSE,
"Elastic Net" = result$en_all[[sh]]$final$CV$CV_RMSE,
GAM = result$gam_all[[sh]]$final$CV$CV_RMSE
)
model_best_aic <- names(aicc_vals)[which.min(aicc_vals)]
model_best_cv <- names(cv_vals)[which.min(cv_vals)]
rows[[sh]] <- data.frame(
Sheet = sh,
AICc_Polynomial = aicc_vals["Polynomial"],
AICc_ElasticNet = aicc_vals["Elastic Net"],
AICc_GAM = aicc_vals["GAM"],
CV_RMSE_Polynomial = cv_vals["Polynomial"],
CV_RMSE_ElasticNet = cv_vals["Elastic Net"],
CV_RMSE_GAM = cv_vals["GAM"],
Model_Terbaik_AIC = model_best_aic,
Model_Terbaik_CV = model_best_cv,
Sama = (model_best_aic == model_best_cv)
)
}
do.call(rbind, rows)
}
best_model_df <- build_best_model_df(result)
agreement_rate <- mean(best_model_df$Sama)
cat(sprintf("Tingkat kesepakatan model terbaik (AICc vs CV) antar sheet: %.2f%%\n", agreement_rate * 100))## Tingkat kesepakatan model terbaik (AICc vs CV) antar sheet: 100.00%best_model_df## Sheet AICc_Polynomial AICc_ElasticNet AICc_GAM CV_RMSE_Polynomial
## Sheet1 Sheet1 22252.78 23318.06 21837.24 4.279636
## Sheet2 Sheet2 22241.77 23265.64 21818.80 4.273885
## Sheet3 Sheet3 22249.39 23304.35 21848.26 4.276527
## Sheet4 Sheet4 22069.54 23250.72 21655.34 4.227976
## Sheet5 Sheet5 22203.59 23306.63 21783.32 4.264502
## CV_RMSE_ElasticNet CV_RMSE_GAM Model_Terbaik_AIC Model_Terbaik_CV Sama
## Sheet1 4.574881 4.151538 GAM GAM TRUE
## Sheet2 4.565893 4.153729 GAM GAM TRUE
## Sheet3 4.577456 4.163720 GAM GAM TRUE
## Sheet4 4.567286 4.117419 GAM GAM TRUE
## Sheet5 4.577479 4.142578 GAM GAM TRUE4.3. Perbandingan dan Pemilihan Model
Perbandingan model dilakukan berdasarkan nilai AICc dan CV_RMSE yang diperoleh dari proses pemodelan. AICc digunakan untuk mempertimbangkan kesesuaian model dan kompleksitasnya, sedangkan CV_RMSE digunakan untuk mengukur kesalahan prediksi melalui 10-fold Cross-Validation.
Hasil pemodelan menunjukkan bahwa GAM dipilih berdasarkan AICc maupun CV_RMSE pada seluruh sheet. Polynomial Regression memiliki nilai CV_RMSE yang lebih tinggi dibandingkan GAM, sedangkan Elastic Net menghasilkan nilai CV_RMSE yang paling tinggi. Dengan demikian, kedua kriteria pemilihan model memberikan hasil yang konsisten dalam menentukan model yang digunakan.
Kesesuaian peringkat model berdasarkan AICc dan CV_RMSE selanjutnya diuji menggunakan korelasi peringkat Spearman. Hasil pengujian menunjukkan koefisien Spearman sebesar 1 pada seluruh sheet dengan p-value sebesar 0,3333. Nilai koefisien tersebut menunjukkan bahwa urutan ketiga model berdasarkan AICc dan CV_RMSE sama pada setiap sheet. Namun, p-value lebih besar dari 0,05 sehingga hubungan tersebut tidak menunjukkan signifikansi statistik pada taraf 5%.
Secara keseluruhan, hasil perbandingan menunjukkan bahwa AICc dan 10-fold Cross-Validation memberikan hasil pemilihan model yang konsisten pada seluruh sheet. Hasil tersebut selanjutnya digunakan sebagai dasar dalam evaluasi model pada data uji.
Syntax yang digunakan:
build_rank_corr_df <- function(best_model_df) {
rows <- list()
for (i in 1:nrow(best_model_df)) {
row <- best_model_df[i, ]
aicc <- c(row$AICc_Polynomial, row$AICc_ElasticNet, row$AICc_GAM)
cv <- c(row$CV_RMSE_Polynomial, row$CV_RMSE_ElasticNet, row$CV_RMSE_GAM)
res <- cor.test(aicc, cv, method = "spearman")
rows[[i]] <- data.frame(
Sheet = row$Sheet,
Spearman_rho = unname(res$estimate),
p_value = res$p.value
)
}
do.call(rbind, rows)
}
rank_corr_df <- build_rank_corr_df(best_model_df)
rank_corr_df## Sheet Spearman_rho p_value
## 1 Sheet1 1 0.3333333
## 2 Sheet2 1 0.3333333
## 3 Sheet3 1 0.3333333
## 4 Sheet4 1 0.3333333
## 5 Sheet5 1 0.33333334.4. Evaluasi Model pada Data Uji
Evaluasi model pada data uji dilakukan untuk mengetahui kemampuan model dalam menghasilkan prediksi pada observasi yang tidak digunakan dalam proses pelatihan. Kinerja model dievaluasi menggunakan RMSE dan R².
Hasil evaluasi menunjukkan bahwa model yang diperoleh mampu memberikan nilai kesalahan prediksi yang relatif rendah pada data uji. Nilai RMSE berada pada kisaran sekitar 4,0, sedangkan nilai R² berada pada kisaran sekitar 0,94. Hasil tersebut menunjukkan bahwa model mampu menjelaskan sebagian besar variasi PE pada data uji.
Evaluasi data uji juga digunakan untuk melihat kesesuaian hasil dengan pemilihan model berdasarkan AICc dan CV_RMSE. Hasil ini menjadi dasar untuk menilai kinerja model sebelum dilakukan analisis perbedaan kinerja antar-sheet.
Syntax yang digunakan:
detail_rows <- list()
for (sh in result$sheet_names) {
for (model_name in c("Polynomial", "Elastic Net", "GAM")) {
all_dict <- switch(model_name,
"Polynomial" = result$poly_all,
"Elastic Net" = result$en_all,
"GAM" = result$gam_all)
for (route in c("AIC", "CV")) {
f <- all_dict[[sh]]$final[[route]]
detail_rows[[length(detail_rows) + 1]] <- data.frame(
Sheet = sh, Model = model_name, Rute = route,
Konfigurasi = f$config, AICc = f$AICc, CV_RMSE = f$CV_RMSE,
Test_RMSE = f$test_RMSE, Test_R2 = f$test_R2
)
}
}
}
detail_df <- do.call(rbind, detail_rows)
detail_df## Sheet Model Rute Konfigurasi AICc CV_RMSE Test_RMSE
## 1 Sheet1 Polynomial AIC degree=3 22252.78 NA 4.116409
## 2 Sheet1 Polynomial CV degree=3 NA 4.279636 4.116409
## 3 Sheet1 Elastic Net AIC lambda=0.0038, l1=0.1 23318.06 NA 4.449722
## 4 Sheet1 Elastic Net CV lambda=0.0038, l1=0.1 NA 4.574881 4.449722
## 5 Sheet1 GAM AIC lambda(GCV)~=0.05372 21837.24 NA 3.992165
## 6 Sheet1 GAM CV lambda=0.004642 NA 4.151538 3.998908
## 7 Sheet2 Polynomial AIC degree=3 22241.77 NA 4.129567
## 8 Sheet2 Polynomial CV degree=3 NA 4.273885 4.129567
## 9 Sheet2 Elastic Net AIC lambda=0.0038, l1=0.1 23265.64 NA 4.512531
## 10 Sheet2 Elastic Net CV lambda=0.0038, l1=0.1 NA 4.565893 4.512531
## 11 Sheet2 GAM AIC lambda(GCV)~=0.05199 21818.80 NA 4.017704
## 12 Sheet2 GAM CV lambda=0.004642 NA 4.153729 4.020662
## 13 Sheet3 Polynomial AIC degree=3 22249.39 NA 4.120109
## 14 Sheet3 Polynomial CV degree=3 NA 4.276527 4.120109
## 15 Sheet3 Elastic Net AIC lambda=0.0038, l1=0.1 23304.35 NA 4.465787
## 16 Sheet3 Elastic Net CV lambda=0.0038, l1=0.1 NA 4.577456 4.465787
## 17 Sheet3 GAM AIC lambda(GCV)~=0.1037 21848.26 NA 3.983542
## 18 Sheet3 GAM CV lambda=0.004642 NA 4.163720 3.983234
## 19 Sheet4 Polynomial AIC degree=3 22069.54 NA 4.322818
## 20 Sheet4 Polynomial CV degree=3 NA 4.227976 4.322818
## 21 Sheet4 Elastic Net AIC lambda=0.0038, l1=0.1 23250.72 NA 4.531404
## 22 Sheet4 Elastic Net CV lambda=0.0038, l1=0.1 NA 4.567286 4.531404
## 23 Sheet4 GAM AIC lambda(GCV)~=0.047 21655.34 NA 4.195046
## 24 Sheet4 GAM CV lambda=0.02154 NA 4.117419 4.208255
## 25 Sheet5 Polynomial AIC degree=3 22203.59 NA 4.171133
## 26 Sheet5 Polynomial CV degree=3 NA 4.264502 4.171133
## 27 Sheet5 Elastic Net AIC lambda=0.0038, l1=0.1 23306.63 NA 4.462542
## 28 Sheet5 Elastic Net CV lambda=0.0038, l1=0.1 NA 4.577479 4.462542
## 29 Sheet5 GAM AIC lambda(GCV)~=0.04472 21783.32 NA 4.050962
## 30 Sheet5 GAM CV lambda=0.004642 NA 4.142578 4.059058
## Test_R2
## 1 0.9414030
## 2 0.9414030
## 3 0.9315294
## 4 0.9315294
## 5 0.9448868
## 6 0.9447005
## 7 0.9433802
## 8 0.9433802
## 9 0.9323917
## 10 0.9323917
## 11 0.9464061
## 12 0.9463272
## 13 0.9419333
## 14 0.9419333
## 15 0.9317809
## 16 0.9317809
## 17 0.9457189
## 18 0.9457273
## 19 0.9340119
## 20 0.9340119
## 21 0.9274901
## 22 0.9274901
## 23 0.9378551
## 24 0.9374632
## 25 0.9404738
## 26 0.9404738
## 27 0.9318659
## 28 0.9318659
## 29 0.9438543
## 30 0.94362964.5. Perbedaan Kinerja Antar Randomisasi
Perbedaan kinerja antar-randomisasi dianalisis menggunakan ANOVA satu arah dan uji Kruskal-Wallis berdasarkan nilai CV-fold RMSE. Kedua pengujian digunakan untuk mengetahui apakah terdapat perbedaan kinerja model antar-sheet.
Hasil ANOVA menunjukkan p-value sebesar 0,9961 untuk Polynomial Regression, 0,99995 untuk Elastic Net, dan 0,9967 untuk GAM. Sementara itu, uji Kruskal-Wallis menghasilkan p-value sebesar 0,8920, 0,9328, dan 0,9067 untuk ketiga metode secara berturut-turut.
Seluruh p-value lebih besar dari 0,05. Dengan demikian, tidak terdapat bukti statistik yang menunjukkan adanya perbedaan kinerja antar-randomisasi berdasarkan nilai CV-fold RMSE, baik menurut ANOVA maupun Kruskal-Wallis.
Hasil ini diperkuat oleh visualisasi boxplot CV-fold RMSE yang menunjukkan pola sebaran yang relatif serupa antar-sheet. Dengan demikian, variasi sheet tidak menghasilkan perbedaan kinerja model yang berarti pada proses Cross-Validation.
cv_fold_store <- list(Polynomial = list(), "Elastic Net" = list(), GAM = list())
for (sh in result$sheet_names) {
deg <- get_config_degree(result$poly_all[[sh]]$final$CV$config)
cv_fold_store$Polynomial[[sh]] <- result$poly_all[[sh]]$cv_fold_rmse[[as.character(deg)]]
cv_fold_store$`Elastic Net`[[sh]] <- result$en_all[[sh]]$cv_fold_rmse
cv_fold_store$GAM[[sh]] <- result$gam_all[[sh]]$cv_fold_rmse
}
anova_rows <- list()
for (model_name in names(cv_fold_store)) {
per_sheet <- cv_fold_store[[model_name]]
df_anova <- do.call(rbind, lapply(names(per_sheet), function(sh) {
data.frame(Sheet = sh, RMSE = per_sheet[[sh]])
}))
aov_res <- summary(aov(RMSE ~ Sheet, data = df_anova))[[1]]
f_stat <- aov_res["Sheet", "F value"]
p_anova <- aov_res["Sheet", "Pr(>F)"]
kw_res <- kruskal.test(RMSE ~ Sheet, data = df_anova)
h_stat <- unname(kw_res$statistic)
p_kw <- kw_res$p.value
anova_rows[[model_name]] <- data.frame(
Model = model_name,
ANOVA_F = f_stat, ANOVA_p = p_anova,
KruskalWallis_H = h_stat, KruskalWallis_p = p_kw,
Signifikan_p_lt_0_05 = (p_anova < 0.05)
)
}
anova_df <- do.call(rbind, anova_rows)
anova_df## Model ANOVA_F ANOVA_p KruskalWallis_H KruskalWallis_p
## Polynomial Polynomial 0.044540876 0.9961140 1.1143529 0.8919883
## Elastic Net Elastic Net 0.004830853 0.9999516 0.8414118 0.9328116
## GAM GAM 0.040660400 0.9967435 1.0202353 0.9067123
## Signifikan_p_lt_0_05
## Polynomial FALSE
## Elastic Net FALSE
## GAM FALSEplot_data_list <- list()
for (model_name in names(cv_fold_store)) {
per_sheet <- cv_fold_store[[model_name]]
for (sh in names(per_sheet)) {
plot_data_list[[length(plot_data_list) + 1]] <- data.frame(
Model = model_name, Sheet = sh, CV_Fold_RMSE = per_sheet[[sh]]
)
}
}
plot_df <- do.call(rbind, plot_data_list)
ggplot(plot_df, aes(x = Sheet, y = CV_Fold_RMSE)) +
geom_boxplot(fill = "lightsteelblue") +
facet_wrap(~ Model, scales = "free_y") +
labs(title = "Sebaran CV-Fold RMSE per Sheet (dasar uji ANOVA/Kruskal-Wallis)") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
4.6 Pembahasan
Hasil eksplorasi menunjukkan bahwa data memiliki hubungan yang kuat antara beberapa variabel prediktor dengan PE, terutama AT dan V. Selain itu, ditemukan adanya indikasi multikolinearitas serta penyimpangan normalitas dan heteroskedastisitas pada residual model linear awal. Kondisi tersebut menunjukkan bahwa hubungan antara prediktor dan PE tidak sepenuhnya sederhana dan linear, sehingga penggunaan metode yang mampu menangkap pola nonlinier menjadi relevan.
Hasil pemodelan menunjukkan bahwa GAM menghasilkan nilai CV_RMSE yang lebih rendah dibandingkan Polynomial Regression dan Elastic Net. Hasil tersebut menunjukkan bahwa pendekatan fungsi smooth pada GAM mampu menangkap pola hubungan nonlinier antara prediktor dan PE dengan lebih baik pada data penelitian. Polynomial Regression juga mampu memodelkan hubungan nonlinier, tetapi bentuk hubungan dibatasi oleh derajat polynomial yang digunakan. Sementara itu, Elastic Net lebih berfokus pada pengendalian kompleksitas dan multikolinearitas melalui penalti sehingga fleksibilitasnya dalam menangkap pola nonlinier lebih terbatas.
Konsistensi hasil pemilihan model berdasarkan AICc dan 10-fold Cross-Validation menunjukkan bahwa kedua kriteria memberikan urutan model yang sama pada setiap randomisasi. Hasil korelasi Spearman yang bernilai 1 juga memperkuat kesamaan peringkat tersebut. Hal ini menunjukkan bahwa keputusan pemilihan model relatif konsisten terhadap kriteria evaluasi yang digunakan.
Hasil pengujian antar-randomisasi menunjukkan tidak terdapat perbedaan kinerja yang signifikan berdasarkan CV-fold RMSE. Kondisi ini dapat dikaitkan dengan kesamaan karakteristik data yang digunakan pada setiap randomisasi. Dengan demikian, perubahan randomisasi tidak menghasilkan perubahan kinerja model yang berarti.
Secara keseluruhan, hasil penelitian menunjukkan bahwa GAM memberikan kemampuan pemodelan yang paling sesuai untuk menangkap hubungan antara AT, V, AP, RH, dan PE pada data penelitian. Hasil evaluasi pada data uji yang menunjukkan nilai RMSE sekitar 4 dan R² sekitar 0,94 juga mendukung kemampuan model dalam menghasilkan prediksi yang baik pada observasi yang tidak digunakan dalam proses pelatihan.
Bab 5. Kesimpulan
Berdasarkan hasil analisis dan pembahasan yang telah dilakukan, diperoleh beberapa kesimpulan sebagai berikut.
Data menunjukkan adanya hubungan yang cukup kuat antara variabel prediktor dengan Net Electrical Output (PE), terutama pada variabel Ambient Temperature (AT) dan Exhaust Vacuum (V). Pemeriksaan awal juga menunjukkan adanya multikolinearitas serta penyimpangan normalitas dan heteroskedastisitas pada residual model linear awal.
Tiga metode yang digunakan, yaitu Polynomial Regression, Elastic Net, dan Generalized Additive Model (GAM), mampu digunakan untuk memodelkan hubungan antara AT, V, AP, RH, dan PE. Berdasarkan hasil 10-fold Cross-Validation, GAM menghasilkan nilai CV_RMSE yang lebih rendah dibandingkan Polynomial Regression dan Elastic Net pada seluruh randomisasi.
Pemilihan model berdasarkan AICc dan CV_RMSE menghasilkan model yang sama pada seluruh randomisasi. Korelasi peringkat Spearman sebesar 1 menunjukkan bahwa urutan model berdasarkan kedua kriteria tersebut identik.
Evaluasi pada data uji menunjukkan bahwa model menghasilkan nilai RMSE sekitar 4 dan R² sekitar 0,94. Hasil tersebut menunjukkan bahwa model mampu menghasilkan prediksi PE dengan tingkat kesalahan yang relatif rendah dan menjelaskan sebagian besar variasi data uji.
Hasil ANOVA dan Kruskal-Wallis menunjukkan tidak terdapat perbedaan kinerja yang signifikan antar-randomisasi berdasarkan nilai CV-fold RMSE. Dengan demikian, perubahan randomisasi pada pembentukan data tidak menunjukkan perubahan kinerja model yang berarti.
Secara keseluruhan, berdasarkan kriteria AICc, 10-fold Cross-Validation, dan evaluasi pada data uji, GAM memberikan hasil pemodelan yang paling sesuai untuk menggambarkan hubungan antara variabel AT, V, AP, RH, dan PE pada penelitian ini.