data_raw <- read.csv(
"C:/Users/HP/OneDrive/berkasku/SEMESTER 5/ANMUL/data set uts anmul.csv",
sep = ";",
header = TRUE,
stringsAsFactors = FALSE
)
cat("DIMENSI DATA AWAL\n")
## DIMENSI DATA AWAL
print(dim(data_raw))
## [1] 649 33
cat("\nNAMA VARIABEL\n")
##
## NAMA VARIABEL
print(names(data_raw))
## [1] "school" "sex" "age" "address" "famsize"
## [6] "Pstatus" "Medu" "Fedu" "Mjob" "Fjob"
## [11] "reason" "guardian" "traveltime" "studytime" "failures"
## [16] "schoolsup" "famsup" "paid" "activities" "nursery"
## [21] "higher" "internet" "romantic" "famrel" "freetime"
## [26] "goout" "Dalc" "Walc" "health" "absences"
## [31] "G1" "G2" "G3"
cat("\nSTRUKTUR DATA\n")
##
## STRUKTUR DATA
str(data_raw)
## 'data.frame': 649 obs. of 33 variables:
## $ school : chr "GP" "GP" "GP" "GP" ...
## $ sex : chr "F" "F" "F" "F" ...
## $ age : int 18 17 15 15 16 16 16 17 15 15 ...
## $ address : chr "U" "U" "U" "U" ...
## $ famsize : chr "GT3" "GT3" "LE3" "GT3" ...
## $ Pstatus : chr "A" "T" "T" "T" ...
## $ Medu : int 4 1 1 4 3 4 2 4 3 3 ...
## $ Fedu : int 4 1 1 2 3 3 2 4 2 4 ...
## $ Mjob : chr "at_home" "at_home" "at_home" "health" ...
## $ Fjob : chr "teacher" "other" "other" "services" ...
## $ reason : chr "course" "course" "other" "home" ...
## $ guardian : chr "mother" "father" "mother" "mother" ...
## $ traveltime: int 2 1 1 1 1 1 1 2 1 1 ...
## $ studytime : int 2 2 2 3 2 2 2 2 2 2 ...
## $ failures : int 0 0 0 0 0 0 0 0 0 0 ...
## $ schoolsup : chr "yes" "no" "yes" "no" ...
## $ famsup : chr "no" "yes" "no" "yes" ...
## $ paid : chr "no" "no" "no" "no" ...
## $ activities: chr "no" "no" "no" "yes" ...
## $ nursery : chr "yes" "no" "yes" "yes" ...
## $ higher : chr "yes" "yes" "yes" "yes" ...
## $ internet : chr "no" "yes" "yes" "yes" ...
## $ romantic : chr "no" "no" "no" "yes" ...
## $ famrel : int 4 5 4 3 4 5 4 4 4 5 ...
## $ freetime : int 3 3 3 2 3 4 4 1 2 5 ...
## $ goout : int 4 3 2 2 2 2 4 4 2 1 ...
## $ Dalc : int 1 1 2 1 1 1 1 1 1 1 ...
## $ Walc : int 1 1 3 1 2 2 1 1 1 1 ...
## $ health : int 3 3 3 5 5 5 3 1 1 5 ...
## $ absences : int 4 2 6 0 0 6 0 2 0 0 ...
## $ G1 : int 0 9 12 14 11 12 13 10 15 12 ...
## $ G2 : int 11 11 13 14 13 12 12 13 16 12 ...
## $ G3 : int 11 11 12 14 13 13 13 13 17 13 ...
cat("\nLIMA BARIS PERTAMA\n")
##
## LIMA BARIS PERTAMA
print(head(data_raw))
## school sex age address famsize Pstatus Medu Fedu Mjob Fjob reason
## 1 GP F 18 U GT3 A 4 4 at_home teacher course
## 2 GP F 17 U GT3 T 1 1 at_home other course
## 3 GP F 15 U LE3 T 1 1 at_home other other
## 4 GP F 15 U GT3 T 4 2 health services home
## 5 GP F 16 U GT3 T 3 3 other other home
## 6 GP M 16 U LE3 T 4 3 services other reputation
## guardian traveltime studytime failures schoolsup famsup paid activities
## 1 mother 2 2 0 yes no no no
## 2 father 1 2 0 no yes no no
## 3 mother 1 2 0 yes no no no
## 4 mother 1 3 0 no yes no yes
## 5 father 1 2 0 no yes no no
## 6 mother 1 2 0 no yes no yes
## nursery higher internet romantic famrel freetime goout Dalc Walc health
## 1 yes yes no no 4 3 4 1 1 3
## 2 no yes yes no 5 3 3 1 1 3
## 3 yes yes yes no 4 3 2 2 3 3
## 4 yes yes yes yes 3 2 2 1 1 5
## 5 yes yes no no 4 3 2 1 2 5
## 6 yes yes yes no 5 4 2 1 2 5
## absences G1 G2 G3
## 1 4 0 11 11
## 2 2 9 11 11
## 3 6 12 13 12
## 4 0 14 14 14
## 5 0 11 13 13
## 6 6 12 12 13
cat("\nJUMLAH MISSING VALUE PER VARIABEL\n")
##
## JUMLAH MISSING VALUE PER VARIABEL
print(colSums(is.na(data_raw)))
## school sex age address famsize Pstatus Medu
## 0 0 0 0 0 0 0
## Fedu Mjob Fjob reason guardian traveltime studytime
## 0 0 0 0 0 0 0
## failures schoolsup famsup paid activities nursery higher
## 0 0 0 0 0 0 0
## internet romantic famrel freetime goout Dalc Walc
## 0 0 0 0 0 0 0
## health absences G1 G2 G3
## 0 0 0 0 0
cat("\nJUMLAH BARIS DUPLIKAT PADA DATA AWAL\n")
##
## JUMLAH BARIS DUPLIKAT PADA DATA AWAL
print(sum(duplicated(data_raw)))
## [1] 0
# ============================================================
# CLEANING DAN PREPROCESSING
# ============================================================
variabel_pilih <- c(
"G3", "studytime", "failures",
"absences", "Medu", "Fedu"
)
# Memastikan semua variabel yang dibutuhkan tersedia
if (!all(variabel_pilih %in% names(data_raw))) {
stop("Ada variabel yang tidak ditemukan. Periksa nama kolom data.")
}
data_pilih <- data_raw[, variabel_pilih, drop = FALSE]
# Konversi ke numerik
data_pilih[] <- lapply(data_pilih, function(x) {
as.numeric(as.character(x))
})
# Hitung data yang hilang setelah konversi
cat("\nMISSING VALUE PADA VARIABEL PENELITIAN\n")
##
## MISSING VALUE PADA VARIABEL PENELITIAN
print(colSums(is.na(data_pilih)))
## G3 studytime failures absences Medu Fedu
## 0 0 0 0 0 0
# Hapus baris dengan missing value
data_bersih <- data_pilih[complete.cases(data_pilih), ]
# Memeriksa rentang nilai yang valid
baris_valid <- with(
data_bersih,
G3 >= 0 & G3 <= 20 &
studytime %in% 1:4 &
failures >= 0 & failures <= 4 &
absences >= 0 &
Medu %in% 0:4 &
Fedu %in% 0:4
)
data_bersih <- data_bersih[baris_valid, , drop = FALSE]
if (nrow(data_bersih) < 10) {
stop("Data valid terlalu sedikit. Periksa format dan isi CSV.")
}
cat("\nJUMLAH BARIS DATA AWAL:", nrow(data_raw), "\n")
##
## JUMLAH BARIS DATA AWAL: 649
cat("JUMLAH BARIS DATA BERSIH:", nrow(data_bersih), "\n")
## JUMLAH BARIS DATA BERSIH: 649
cat("JUMLAH BARIS YANG TIDAK DIGUNAKAN:",
nrow(data_raw) - nrow(data_bersih), "\n")
## JUMLAH BARIS YANG TIDAK DIGUNAKAN: 0
cat("\nDUPLIKASI PADA VARIABEL TERPILIH:",
sum(duplicated(data_bersih)), "\n")
##
## DUPLIKASI PADA VARIABEL TERPILIH: 103
print(summary(data_bersih))
## G3 studytime failures absences
## Min. : 0.00 Min. :1.000 Min. :0.0000 Min. : 0.000
## 1st Qu.:10.00 1st Qu.:1.000 1st Qu.:0.0000 1st Qu.: 0.000
## Median :12.00 Median :2.000 Median :0.0000 Median : 2.000
## Mean :11.91 Mean :1.931 Mean :0.2219 Mean : 3.659
## 3rd Qu.:14.00 3rd Qu.:2.000 3rd Qu.:0.0000 3rd Qu.: 6.000
## Max. :19.00 Max. :4.000 Max. :3.0000 Max. :32.000
## Medu Fedu
## Min. :0.000 Min. :0.000
## 1st Qu.:2.000 1st Qu.:1.000
## Median :2.000 Median :2.000
## Mean :2.515 Mean :2.307
## 3rd Qu.:4.000 3rd Qu.:3.000
## Max. :4.000 Max. :4.000
write.csv(
data_bersih,
"output/student_por_bersih.csv",
row.names = FALSE
)
# ============================================================
# EXPLORATORY DATA ANALYSIS (EDA)
# ============================================================
statistik_deskriptif <- data.frame(
Variabel = names(data_bersih),
Jumlah = sapply(data_bersih, length),
Mean = sapply(data_bersih, mean),
Median = sapply(data_bersih, median),
Standar_Deviasi = sapply(data_bersih, sd),
Minimum = sapply(data_bersih, min),
Maksimum = sapply(data_bersih, max)
)
cat("\nSTATISTIK DESKRIPTIF\n")
##
## STATISTIK DESKRIPTIF
print(statistik_deskriptif, row.names = FALSE)
## Variabel Jumlah Mean Median Standar_Deviasi Minimum Maksimum
## G3 649 11.9060092 12 3.2306562 0 19
## studytime 649 1.9306626 2 0.8295096 1 4
## failures 649 0.2218798 0 0.5932351 0 3
## absences 649 3.6594761 2 4.6407588 0 32
## Medu 649 2.5146379 2 1.1345520 0 4
## Fedu 649 2.3066256 2 1.0999309 0 4
# Histogram nilai akhir
hist(
data_bersih$G3,
main = "Distribusi Nilai Akhir Bahasa Portugis",
xlab = "Nilai Akhir (G3)",
ylab = "Frekuensi",
col = "#F4B6C2",
border = "white",
breaks = seq(-0.5, 20.5, by = 1)
)
# Boxplot nilai akhir
boxplot(
data_bersih$G3,
main = "Sebaran Nilai Akhir Siswa",
ylab = "Nilai Akhir (G3)",
col = "#C9B6E4",
border = "#8E78B7"
)
# Nilai akhir berdasarkan waktu belajar
boxplot(
G3 ~ factor(studytime),
data = data_bersih,
main = "Nilai Akhir Berdasarkan Waktu Belajar",
xlab = "Waktu Belajar per Minggu",
ylab = "Nilai Akhir (G3)",
names = c("<2 jam", "2-5 jam", "5-10 jam", ">10 jam"),
col = c("#F4B6C2", "#F8D7A3", "#B9DCC4", "#B8D4F0")
)
# Nilai akhir berdasarkan jumlah kegagalan kelas sebelumnya
boxplot(
G3 ~ factor(failures),
data = data_bersih,
main = "Nilai Akhir Berdasarkan Kegagalan Kelas",
xlab = "Jumlah Kegagalan Sebelumnya",
ylab = "Nilai Akhir (G3)",
col = "#D6C6F0"
)
# Ketidakhadiran dan nilai akhir
plot(
data_bersih$absences,
data_bersih$G3,
main = "Ketidakhadiran dan Nilai Akhir",
xlab = "Jumlah Ketidakhadiran",
ylab = "Nilai Akhir (G3)",
pch = 19,
col = "#91B8E8"
)
# Matriks korelasi untuk eksplorasi awal
cat("\nMATRIKS KORELASI\n")
##
## MATRIKS KORELASI
print(round(cor(data_bersih), 3))
## G3 studytime failures absences Medu Fedu
## G3 1.000 0.250 -0.393 -0.091 0.240 0.212
## studytime 0.250 1.000 -0.147 -0.118 0.097 0.050
## failures -0.393 -0.147 1.000 0.123 -0.172 -0.166
## absences -0.091 -0.118 0.123 1.000 -0.009 0.030
## Medu 0.240 0.097 -0.172 -0.009 1.000 0.647
## Fedu 0.212 0.050 -0.166 0.030 0.647 1.000
# ============================================================
# REGRESI LINEAR BERGANDA (MODEL AWAL)
# ============================================================
model_mlr <- lm(
G3 ~ studytime + failures + absences + Medu + Fedu,
data = data_bersih
)
cat("\nRINGKASAN MODEL MLR AWAL\n")
##
## RINGKASAN MODEL MLR AWAL
print(summary(model_mlr))
##
## Call:
## lm(formula = G3 ~ studytime + failures + absences + Medu + Fedu,
## data = data_bersih)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.2913 -1.4828 0.0462 1.7584 6.7766
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 9.67420 0.41028 23.580 < 2e-16 ***
## studytime 0.71017 0.13811 5.142 3.61e-07 ***
## failures -1.79892 0.19584 -9.186 < 2e-16 ***
## absences -0.02113 0.02457 -0.860 0.3900
## Medu 0.33845 0.13083 2.587 0.0099 **
## Fedu 0.21074 0.13466 1.565 0.1181
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.861 on 643 degrees of freedom
## Multiple R-squared: 0.222, Adjusted R-squared: 0.2159
## F-statistic: 36.69 on 5 and 643 DF, p-value: < 2.2e-16
# ============================================================
# SELEKSI VARIABEL DENGAN AIC
# ============================================================
model_aic <- step(
model_mlr,
direction = "both",
trace = 0
)
cat("\nMODEL TERPILIH BERDASARKAN AIC\n")
##
## MODEL TERPILIH BERDASARKAN AIC
print(formula(model_aic))
## G3 ~ studytime + failures + Medu + Fedu
cat("\nRINGKASAN MODEL AIC\n")
##
## RINGKASAN MODEL AIC
print(summary(model_aic))
##
## Call:
## lm(formula = G3 ~ studytime + failures + Medu + Fedu, data = data_bersih)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.2116 -1.5352 0.0557 1.7884 6.8962
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 9.5868 0.3974 24.123 < 2e-16 ***
## studytime 0.7223 0.1374 5.258 1.98e-07 ***
## failures -1.8179 0.1946 -9.344 < 2e-16 ***
## Medu 0.3405 0.1308 2.604 0.00944 **
## Fedu 0.2045 0.1344 1.521 0.12863
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.86 on 644 degrees of freedom
## Multiple R-squared: 0.2211, Adjusted R-squared: 0.2163
## F-statistic: 45.7 on 4 and 644 DF, p-value: < 2.2e-16
cat("\nPERBANDINGAN AIC MODEL\n")
##
## PERBANDINGAN AIC MODEL
print(AIC(model_mlr, model_aic))
## df AIC
## model_mlr 7 3214.021
## model_aic 6 3212.767
# ============================================================
# UJI ASUMSI MODEL TERPILIH
# ============================================================
if (!requireNamespace("car", quietly = TRUE)) {
}
if (!requireNamespace("lmtest", quietly = TRUE)) {
}
if (!requireNamespace("sandwich", quietly = TRUE)) {
}
# Normalitas residual
cat("\nUJI NORMALITAS RESIDUAL\n")
##
## UJI NORMALITAS RESIDUAL
if (length(residuals(model_aic)) <= 5000) {
print(shapiro.test(residuals(model_aic)))
} else {
cat("Jumlah residual melebihi batas Shapiro-Wilk.\n")
}
##
## Shapiro-Wilk normality test
##
## data: residuals(model_aic)
## W = 0.93736, p-value = 7.048e-16
# Multikolinearitas
cat("\nUJI MULTIKOLINEARITAS (VIF)\n")
##
## UJI MULTIKOLINEARITAS (VIF)
print(car::vif(model_aic))
## studytime failures Medu Fedu
## 1.028488 1.055210 1.744130 1.732135
# Heteroskedastisitas
cat("\nUJI BREUSCH-PAGAN\n")
##
## UJI BREUSCH-PAGAN
print(lmtest::bptest(model_aic))
##
## studentized Breusch-Pagan test
##
## data: model_aic
## BP = 13.032, df = 4, p-value = 0.01112
# ============================================================
# ROBUST STANDARD ERRORS
# ============================================================
cat("\nKOEFISIEN DENGAN ROBUST STANDARD ERRORS HC3\n")
##
## KOEFISIEN DENGAN ROBUST STANDARD ERRORS HC3
hasil_robust <- lmtest::coeftest(
model_aic,
vcov. = sandwich::vcovHC(model_aic, type = "HC3")
)
print(hasil_robust)
##
## t test of coefficients:
##
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 9.58679 0.38977 24.5962 < 2.2e-16 ***
## studytime 0.72230 0.13301 5.4305 7.971e-08 ***
## failures -1.81795 0.22123 -8.2176 1.142e-15 ***
## Medu 0.34051 0.13826 2.4628 0.01405 *
## Fedu 0.20454 0.14244 1.4360 0.15150
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# ============================================================
# GRAFIK DIAGNOSTIK MODEL TERPILIH
# ============================================================
par(
mfrow = c(2, 2),
bg = "white",
mar = c(4, 4, 3, 1)
)
plot(model_aic, which = 1, col = "#E78FB3", pch = 19)
plot(model_aic, which = 2, col = "#91B8E8", pch = 19)
plot(model_aic, which = 3, col = "#A8CFA8", pch = 19)
plot(model_aic, which = 5, col = "#C3A6E8", pch = 19)
par(mfrow = c(1, 1))