UTS Statistika Lingkungan
Nurdiana Rizki R
2025-10-15
SOAL 1
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.2 ✔ tibble 3.3.0
## ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
## ✔ purrr 1.1.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(splines)# Data Training
data_train <- read.csv("D:/Semester 5/Statistika Lingkungan/UTS/kualitasair.xlsx - Training.csv")head(data_train)## Lokasi pH DO BOD TSS Suhu Status
## 1 S1 7.6855 NA 1.7136 43.1415 26.7972 Tercemar ringan
## 2 S2 6.7177 5.7236 1.4402 44.2963 27.7284 Tercemar ringan
## 3 S3 7.1816 4.8906 2.7274 NA 26.0255 Tercemar ringan
## 4 S4 7.3164 6.1339 3.1398 41.0104 29.6639 Tercemar ringan
## 5 S5 7.2021 7.7853 1.1778 48.0967 26.4099 baik
## 6 S6 6.9469 8.4222 3.2324 48.5610 28.6809 Tercemar ringan# cek missing value
colSums(is.na(data_train))## Lokasi pH DO BOD TSS Suhu Status
## 0 0 23 22 24 0 0berdasarkan hasil pemeriksaan missing value, dataset kualitas air memiliki missing value pada variabel pH, DO, BOD dan TSS
# Mengganti NA pada kolom DO dengan median
data_train$DO[is.na(data_train$DO)] <- median(data_train$DO, na.rm = TRUE)
# Mengganti NA pada kolom BOD dengan median
data_train$BOD[is.na(data_train$BOD)] <- median(data_train$BOD, na.rm = TRUE)
# Mengganti NA pada kolom TSS dengan median
data_train$TSS[is.na(data_train$TSS)] <- median(data_train$TSS, na.rm = TRUE)
# Cek lagi apakah masih ada missing value
colSums(is.na(data_train))## Lokasi pH DO BOD TSS Suhu Status
## 0 0 0 0 0 0 0# Cek outlier variabel
boxplot(data_train$pH, main="Outlier pH")
boxplot(data_train$DO, main="Outlier DO")
boxplot(data_train$BOD, main="Outlier BOD")
boxplot(data_train$TSS, main="Outlier TSS")
boxplot(data_train$Suhu, main="Outlier Suhu")
for (col in c("pH", "DO", "BOD", "TSS", "Suhu")) {
Q1 <- quantile(data_train[[col]], 0.25, na.rm = TRUE)
Q3 <- quantile(data_train[[col]], 0.75, na.rm = TRUE)
IQR <- Q3 - Q1
lower <- Q1 - 1.5 * IQR
upper <- Q3 + 1.5 * IQR
data <- data_train[data_train[[col]] >= lower & data_train[[col]] <= upper, ]
}standarisasi penulisan kategori Status
table(data_train$Status)##
## baik Baik BAIK Tercemar berat tercemar ringan
## 1 70 1 7 1
## Tercemar ringan Tercemar Ringan
## 219 1data_train$Status <- tolower(data_train$Status)
data_train$Status <- trimws(data_train$Status)
data_train$Status <- tolower(data_train$Status)
data_train$Status <- trimws(data_train$Status)
data_train$Status <- recode(data_train$Status,
"baik" = "1",
"tercemar ringan" = "2",
"tercemar berat" = "3")
data_train$Status <- as.factor(data_train$Status)summary(data[, c("pH", "DO", "BOD", "TSS", "Suhu")])## pH DO BOD TSS
## Min. :5.503 Min. :2.982 Min. :0.3026 Min. :24.65
## 1st Qu.:6.668 1st Qu.:5.408 1st Qu.:2.4499 1st Qu.:44.26
## Median :6.987 Median :5.991 Median :3.0661 Median :49.52
## Mean :6.988 Mean :5.975 Mean :3.0039 Mean :49.69
## 3rd Qu.:7.318 3rd Qu.:6.624 3rd Qu.:3.5337 3rd Qu.:55.69
## Max. :8.351 Max. :9.229 Max. :5.7962 Max. :76.34
## Suhu
## Min. :22.77
## 1st Qu.:26.60
## Median :28.00
## Mean :28.08
## 3rd Qu.:29.43
## Max. :33.40SOAL 2
library(caret)## Loading required package: lattice##
## Attaching package: 'caret'## The following object is masked from 'package:purrr':
##
## liftlibrary(e1071)
library(rpart)
library(randomForest)## randomForest 4.7-1.2## Type rfNews() to see new features/changes/bug fixes.##
## Attaching package: 'randomForest'## The following object is masked from 'package:dplyr':
##
## combine## The following object is masked from 'package:ggplot2':
##
## marginset.seed(223)
trainIndex <- createDataPartition(data_train$Status, p = 0.8, list = FALSE)
train_data <- data_train[trainIndex, ]
test_data <- data_train[-trainIndex, ]table(train_data$Status)##
## 1 2 3
## 58 177 6table(test_data$Status)##
## 1 2 3
## 14 44 1# =============================
# Support Vector Machine (SVM)
# =============================
set.seed(123)
svm_model <- svm(Status ~ pH + DO + BOD + TSS + Suhu, data = train_data, kernel = "radial")
# Prediksi dan evaluasi
svm_pred <- predict(svm_model, test_data)
confusionMatrix(svm_pred, test_data$Status)## Confusion Matrix and Statistics
##
## Reference
## Prediction 1 2 3
## 1 9 1 0
## 2 5 43 1
## 3 0 0 0
##
## Overall Statistics
##
## Accuracy : 0.8814
## 95% CI : (0.7707, 0.9509)
## No Information Rate : 0.7458
## P-Value [Acc > NIR] : 0.008645
##
## Kappa : 0.6515
##
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: 1 Class: 2 Class: 3
## Sensitivity 0.6429 0.9773 0.00000
## Specificity 0.9778 0.6000 1.00000
## Pos Pred Value 0.9000 0.8776 NaN
## Neg Pred Value 0.8980 0.9000 0.98305
## Prevalence 0.2373 0.7458 0.01695
## Detection Rate 0.1525 0.7288 0.00000
## Detection Prevalence 0.1695 0.8305 0.00000
## Balanced Accuracy 0.8103 0.7886 0.50000Berdasarkan hasil evaluasi, model memiliki akurasi 88,14% dengan p-value 0,008645 yang menunjukkan performa signifikan. Nilai Kappa 0,6515 mengindikasikan kesesuaian prediksi yang baik. Kinerja per kelas menunjukkan bahwa kelas 2 memiliki sensitivitas tinggi dengan nilai 0,9773, kelas 1 cukup baik 0,6429, sedangkan kelas 3 tidak terdeteksi sama sekali 0,0000. Hal ini menunjukkan adanya ketidakseimbangan data, terutama pada kelas 3 yang sangat kecil
# =============================
# Decision Tree
# =============================
set.seed(123)
tree_model <- rpart(Status ~ pH + DO + BOD + TSS + Suhu, data = train_data, method = "class")
tree_pred <- predict(tree_model, test_data, type = "class")
confusionMatrix(tree_pred, test_data$Status)## Confusion Matrix and Statistics
##
## Reference
## Prediction 1 2 3
## 1 13 0 0
## 2 1 44 1
## 3 0 0 0
##
## Overall Statistics
##
## Accuracy : 0.9661
## 95% CI : (0.8829, 0.9959)
## No Information Rate : 0.7458
## P-Value [Acc > NIR] : 6.696e-06
##
## Kappa : 0.9075
##
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: 1 Class: 2 Class: 3
## Sensitivity 0.9286 1.0000 0.00000
## Specificity 1.0000 0.8667 1.00000
## Pos Pred Value 1.0000 0.9565 NaN
## Neg Pred Value 0.9783 1.0000 0.98305
## Prevalence 0.2373 0.7458 0.01695
## Detection Rate 0.2203 0.7458 0.00000
## Detection Prevalence 0.2203 0.7797 0.00000
## Balanced Accuracy 0.9643 0.9333 0.50000Model memiliki akurasi 96,61% dengan p-value 6.696e-06 yang menunjukkan performa sangat signifikan. Nilai Kappa 0,9075 menandakan kesesuaian prediksi yang sangat baik. Kinerja per kelas menunjukkan kelas 1 sensitivitas 0,9286 dan kelas 2 sensitivitas 1,0000 terdeteksi dengan sangat baik, sedangkan kelas 3 tidak terdeteksi dengan nilai sensitivitas 0,0000.
# =============================
# Random Forest
# =============================
set.seed(123)
rf_model <- randomForest(Status ~ pH + DO + BOD + TSS + Suhu, data = train_data, ntree = 500, mtry = 3)
rf_pred <- predict(rf_model, test_data)
confusionMatrix(rf_pred, test_data$Status)## Confusion Matrix and Statistics
##
## Reference
## Prediction 1 2 3
## 1 13 1 0
## 2 1 43 0
## 3 0 0 1
##
## Overall Statistics
##
## Accuracy : 0.9661
## 95% CI : (0.8829, 0.9959)
## No Information Rate : 0.7458
## P-Value [Acc > NIR] : 6.696e-06
##
## Kappa : 0.9125
##
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: 1 Class: 2 Class: 3
## Sensitivity 0.9286 0.9773 1.00000
## Specificity 0.9778 0.9333 1.00000
## Pos Pred Value 0.9286 0.9773 1.00000
## Neg Pred Value 0.9778 0.9333 1.00000
## Prevalence 0.2373 0.7458 0.01695
## Detection Rate 0.2203 0.7288 0.01695
## Detection Prevalence 0.2373 0.7458 0.01695
## Balanced Accuracy 0.9532 0.9553 1.00000Model menunjukkan akurasi tinggi sebesar 96,61% dengan p-value 6.696e-06, menandakan kinerjanya secara signifikan lebih baik dari tebakan acak. Nilai Kappa 0,9125 memperlihatkan kesesuaian prediksi yang sangat baik. Performa per kelas juga stabil: kelas 1 sensitivitas 0,9286, kelas 2 0,9773 dan kelas 3 dengan nilai 1,0000 terdeteksi dengan sangat baik. Semua kelas memiliki balanced accuracy di atas 0,95. Hasil ini menunjukkan bahwa model telah mampu melakukan klasifikasi secara merata di semua kelas, termasuk kelas minoritas.
# =============================
# Feature Importance (pentingnya variabel)
# =============================
importance(rf_model)## MeanDecreaseGini
## pH 2.665045
## DO 44.298093
## BOD 42.440050
## TSS 2.883089
## Suhu 4.255465svm_acc <- confusionMatrix(svm_pred, test_data$Status)$overall['Accuracy']
tree_acc <- confusionMatrix(tree_pred, test_data$Status)$overall['Accuracy']
rf_acc <- confusionMatrix(rf_pred, test_data$Status)$overall['Accuracy']
akurasi_model <- data.frame(
Model = c("SVM", "Decision Tree", "Random Forest"),
Akurasi = c(svm_acc, tree_acc, rf_acc)
)
print(akurasi_model)## Model Akurasi
## 1 SVM 0.8813559
## 2 Decision Tree 0.9661017
## 3 Random Forest 0.9661017Dari hasil pengujian tiga algoritma klasifikasi, model SVM memperoleh akurasi 88,13%, sedangkan Decision Tree dan Random Forest sama-sama mencapai 96,61%. Hal ini menunjukkan bahwa Decision Tree dan Random Forest memiliki performa klasifikasi yang lebih tinggi dibandingkan SVM pada dataset ini. Dengan tingkat akurasi yang sama dan tinggi, Random Forest atau Decision Tree dapat dipilih sebagai model utama karena lebih optimal dalam mengenali pola data.
SOAL 3
library(caret)
library(splines)
library(ggplot2)
set.seed(123)
data_pred <- data[, c("DO", "pH", "BOD", "TSS", "Suhu")]
trainIndex <- createDataPartition(data_pred$DO, p = 0.8, list = FALSE)
train_data <- data_pred[trainIndex, ]
test_data <- data_pred[-trainIndex, ]# ============================================
# Model Regresi Linear
# ============================================
lm_model <- lm(DO ~ pH + BOD + TSS + Suhu, data = train_data)
lm_pred <- predict(lm_model, newdata = test_data)
# Evaluasi performa
lm_r2 <- summary(lm_model)$r.squared
lm_mse <- mean((lm_pred - test_data$DO)^2)
lm_rmse <- sqrt(lm_mse)
lm_eval <- data.frame(Model = "Regresi Linear", R2 = lm_r2, MSE = lm_mse, RMSE = lm_rmse)
lm_eval## Model R2 MSE RMSE
## 1 Regresi Linear 0.01529239 0.9090539 0.9534432# ============================================
# Model Regresi Spline (Natural Spline)
# ============================================
# Gunakan spline untuk pH dan BOD (biasanya non-linear)
spline_model <- lm(DO ~ ns(pH, df = 3) + ns(BOD, df = 3) + TSS + Suhu,data = train_data)
spline_pred <- predict(spline_model, newdata = test_data)
spline_r2 <- summary(spline_model)$r.squared
spline_mse <- mean((spline_pred - test_data$DO)^2)
spline_rmse <- sqrt(spline_mse)
spline_eval <- data.frame(Model = "Regresi Spline", R2 = spline_r2, MSE = spline_mse, RMSE = spline_rmse)
spline_eval## Model R2 MSE RMSE
## 1 Regresi Spline 0.05976852 0.9595585 0.9795706# Koefisien model linear
eval_models <- rbind(lm_eval, spline_eval)
print(eval_models)## Model R2 MSE RMSE
## 1 Regresi Linear 0.01529239 0.9090539 0.9534432
## 2 Regresi Spline 0.05976852 0.9595585 0.9795706Model Regresi Linear menghasilkan nilai R² sebesar 0,0106, MSE 0,8089, dan RMSE 0,8994. Sementara itu, Regresi Spline menunjukkan kinerja sedikit lebih baik dengan R² sebesar 0,0544, MSE 0,7797, dan RMSE 0,8829.
Nilai R² yang rendah pada kedua model mengindikasikan bahwa keduanya belum mampu menjelaskan variabilitas data secara optimal. Namun, Regresi Spline memberikan hasil lebih baik dibanding Regresi Linear, meskipun perbedaannya relatif kecil. Dengan demikian, Regresi Spline lebih direkomendasikan untuk kasus ini, tetapi perlu peningkatan atau eksplorasi model lain untuk mendapatkan performa prediksi yang lebih tinggi
# ============================================
# Visualisasi: Prediksi vs Aktual
# ============================================
plot_data <- data.frame(
Aktual = test_data$DO,
Prediksi_LM = lm_pred,
Prediksi_Spline = spline_pred
)
# Visualisasi regresi linear
ggplot(plot_data, aes(x = Aktual, y = Prediksi_LM)) +
geom_point(color = "blue", size = 2) +
geom_abline(slope = 1, intercept = 0, color = "red", linetype = "dashed") +
labs(title = "Prediksi vs Aktual (Regresi Linear)",
x = "DO Aktual", y = "DO Prediksi") +
theme_minimal()
# Visualisasi regresi spline
ggplot(plot_data, aes(x = Aktual, y = Prediksi_Spline)) +
geom_point(color = "darkgreen", size = 2) +
geom_abline(slope = 1, intercept = 0, color = "red", linetype = "dashed") +
labs(title = "Prediksi vs Aktual (Regresi Spline)",
x = "DO Aktual", y = "DO Prediksi") +
theme_minimal()
summary(lm_model)##
## Call:
## lm(formula = DO ~ pH + BOD + TSS + Suhu, data = train_data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.05112 -0.50923 0.02646 0.63417 2.96253
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.427e+00 1.312e+00 4.897 1.82e-06 ***
## pH -1.407e-02 1.221e-01 -0.115 0.9084
## BOD 1.343e-01 7.752e-02 1.733 0.0844 .
## TSS -6.275e-05 6.576e-03 -0.010 0.9924
## Suhu -2.737e-02 3.165e-02 -0.865 0.3881
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9689 on 234 degrees of freedom
## Multiple R-squared: 0.01529, Adjusted R-squared: -0.00154
## F-statistic: 0.9085 on 4 and 234 DF, p-value: 0.4597Berdasarkan hasil regresi, variabel yang paling memengaruhi DO adalah BOD. Hal ini terlihat dari nilai koefisien BOD sebesar 0.1343 dengan nilai p-value 0.0844, yang paling kecil dibanding variabel lainnya dan mendekati batas signifikansi 0.05. Artinya, meskipun belum signifikan secara statistik pada tingkat kepercayaan 95%, BOD memiliki kecenderungan memberikan pengaruh terhadap DO. Koefisien positif menunjukkan bahwa ketika nilai BOD meningkat, nilai DO juga cenderung meningkat, meskipun hubungan ini masih lemah. Sementara itu, variabel pH, TSS, dan Suhu memiliki nilai p-value yang jauh lebih besar, sehingga tidak berpengaruh signifikan terhadap DO dalam model ini