# Load library
library(car) # Untuk VIF, Durbin-Watson, dll.
library(lmtest) # Untuk uji Breusch-Pagan, Durbin-Watson
library(nortest) # Untuk uji normalitas (Anderson-Darling)
library(MASS) # Untuk Box-Cox transformation
library(olsrr) # Untuk uji asumsi OLS lengkap
library(ggplot2) # Untuk visualisasi
library(corrplot) # Untuk matriks korelasi
# IMPORT DATA
data <- read.csv("C:/Users/MyBook Hype AMD/Downloads/bodyfat (2).csv")
# PEMISALAN VARIABEL
Y <- data$BodyFat
X1 <- data$Abdomen
X2 <- data$Weight
X3 <- data$Hip
#Buat data frame
data <- data.frame(Y = Y, X1 = X1, X2 = X2, X3 = X3)
# Tampilkan 6 baris pertama
cat("=" %+% rep("=", 60) %+% "\n")
print(head(data))
## Y X1 X2 X3
## 1 12.3 85.2 154.25 94.5
## 2 6.1 83.0 173.25 98.7
## 3 25.3 87.9 154.00 99.2
## 4 10.4 86.4 184.75 101.2
## 5 28.7 100.0 184.25 101.9
## 6 20.9 94.4 210.25 107.8
n<-nrow(data)
cat("banyak data",n) #banyak data
## banyak data 252
# ------------------------------------------------------------
# 2. STATISTIK DESKRIPTIF DAN EKSPLORASI DATA
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
# Statistik deskriptif
summary(data)
## Y X1 X2 X3
## Min. : 0.00 Min. : 69.40 Min. :118.5 Min. : 85.0
## 1st Qu.:12.47 1st Qu.: 84.58 1st Qu.:159.0 1st Qu.: 95.5
## Median :19.20 Median : 90.95 Median :176.5 Median : 99.3
## Mean :19.15 Mean : 92.56 Mean :178.9 Mean : 99.9
## 3rd Qu.:25.30 3rd Qu.: 99.33 3rd Qu.:197.0 3rd Qu.:103.5
## Max. :47.50 Max. :148.10 Max. :363.1 Max. :147.7
# Matriks korelasi
cat("\nMatriks Korelasi:\n")
##
## Matriks Korelasi:
cor_matrix <- cor(data)
print(round(cor_matrix, 3))
## Y X1 X2 X3
## Y 1.000 0.813 0.612 0.625
## X1 0.813 1.000 0.888 0.874
## X2 0.612 0.888 1.000 0.941
## X3 0.625 0.874 0.941 1.000
# Visualisasi matriks korelasi
corrplot(cor_matrix, method = "color", type = "upper",
addCoef.col = "black", tl.col = "black", tl.srt = 45,
title = "Matriks Korelasi Antar Variabel", mar = c(0,0,2,0))
# Scatter plot matrix
pairs(data, main = "Scatter Plot Matrix", pch = 19, col = "steelblue")
# ------------------------------------------------------------
# 3. ESTIMASI MODEL MULTIPLE LINEAR REGRESSION
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
# Fit model MLR
model <- lm(Y ~ X1 + X2 + X3, data = data)
# Ringkasan model
cat("\nRingkasan Model:\n")
##
## Ringkasan Model:
print(summary(model))
##
## Call:
## lm(formula = Y ~ X1 + X2 + X3, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11.6340 -3.2033 -0.0316 3.1601 10.5113
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -45.845537 7.058594 -6.495 4.50e-10 ***
## X1 0.989741 0.058657 16.873 < 2e-16 ***
## X2 -0.147632 0.030866 -4.783 2.97e-06 ***
## X3 -0.001953 0.119859 -0.016 0.987
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.465 on 248 degrees of freedom
## Multiple R-squared: 0.7188, Adjusted R-squared: 0.7154
## F-statistic: 211.3 on 3 and 248 DF, p-value: < 2.2e-16
# Ekstrak koefisien
cat("\nKoefisien Regresi:\n")
##
## Koefisien Regresi:
coef_df <- data.frame(
Variabel = names(coef(model)),
Koefisien = round(coef(model), 4),
Std_Error = round(summary(model)$coefficients[, 2], 4),
t_value = round(summary(model)$coefficients[, 3], 4),
p_value = round(summary(model)$coefficients[, 4], 4)
)
print(coef_df)
## Variabel Koefisien Std_Error t_value p_value
## (Intercept) (Intercept) -45.8455 7.0586 -6.4950 0.000
## X1 X1 0.9897 0.0587 16.8734 0.000
## X2 X2 -0.1476 0.0309 -4.7830 0.000
## X3 X3 -0.0020 0.1199 -0.0163 0.987
# Interval kepercayaan 95%
cat("\nInterval Kepercayaan 95%:\n")
##
## Interval Kepercayaan 95%:
ci <- confint(model, level = 0.95)
print(round(ci, 4))
## 2.5 % 97.5 %
## (Intercept) -59.7480 -31.9431
## X1 0.8742 1.1053
## X2 -0.2084 -0.0868
## X3 -0.2380 0.2341
# ------------------------------------------------------------
# 4. UJI SIGNIFIKANSI SIMULTAN (UJI F)
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
# Ekstrak informasi F-test
f_stat <- summary(model)$fstatistic
f_value <- f_stat[1]
df1 <- f_stat[2]
df2 <- f_stat[3]
f_pvalue <- pf(f_value, df1, df2, lower.tail = FALSE)
cat("\nHipotesis:\n")
##
## Hipotesis:
cat("H0: β1 = β2 = β3 = 0 (Model tidak signifikan)\n")
## H0: β1 = β2 = β3 = 0 (Model tidak signifikan)
cat("H1: Minimal ada satu βj ≠ 0 (Model signifikan)\n")
## H1: Minimal ada satu βj ≠ 0 (Model signifikan)
cat("\nHasil Uji F:\n")
##
## Hasil Uji F:
cat("F-statistic:", round(f_value, 4), "\n")
## F-statistic: 211.3099
cat("df1:", df1, "\n")
## df1: 3
cat("df2:", df2, "\n")
## df2: 248
cat("p-value:", format(f_pvalue, scientific = TRUE), "\n")
## p-value: 5.101077e-68
if (f_pvalue < 0.05) {
cat("Kesimpulan: Tolak H0 → Model signifikan pada α = 5%\n")
} else {
cat("Kesimpulan: Gagal Tolak H0 → Model tidak signifikan\n")
}
## Kesimpulan: Tolak H0 → Model signifikan pada α = 5%
# ------------------------------------------------------------
# 5. UJI SIGNIFIKANSI PARSIAL (UJI T)
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
cat("\nHipotesis untuk setiap βj:\n")
##
## Hipotesis untuk setiap βj:
cat("H0: βj = 0 (Variabel tidak signifikan)\n")
## H0: βj = 0 (Variabel tidak signifikan)
cat("H1: βj ≠ 0 (Variabel signifikan)\n")
## H1: βj ≠ 0 (Variabel signifikan)
cat("\nHasil Uji T:\n")
##
## Hasil Uji T:
for (i in 2:nrow(coef_df)) {
var_name <- coef_df$Variabel[i]
t_val <- coef_df$t_value[i]
p_val <- coef_df$p_value[i]
sig <- ifelse(p_val < 0.001, "Sangat signifikan (***)",
ifelse(p_val < 0.01, "Cukup signifikan (**)",
ifelse(p_val < 0.05, "Signifikan (*)",
ifelse(p_val < 0.1, ".", "Tidak signifikan(ns)"))))
cat(sprintf("%-10s: t = %7.4f, p = %8.4f %s\n", var_name, t_val, p_val, sig))
}
## X1 : t = 16.8734, p = 0.0000 Sangat signifikan (***)
## X2 : t = -4.7830, p = 0.0000 Sangat signifikan (***)
## X3 : t = -0.0163, p = 0.9870 Tidak signifikan(ns)
# ------------------------------------------------------------
# 6. KOEFISIEN DETERMINASI (R² dan Adjusted R²)
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
r_squared <- summary(model)$r.squared
adj_r_squared <- summary(model)$adj.r.squared
resid_se <- summary(model)$sigma
cat("\nR-squared:", round(r_squared, 4), "\n")
##
## R-squared: 0.7188
cat("Adjusted R-squared:", round(adj_r_squared, 4), "\n")
## Adjusted R-squared: 0.7154
cat("Residual Standard Error:", round(resid_se, 4), "\n")
## Residual Standard Error: 4.4646
cat("Interpretasi: Model mampu menjelaskan",
round(r_squared * 100, 2), "% variasi pada Y\n")
## Interpretasi: Model mampu menjelaskan 71.88 % variasi pada Y
cat("\n============================================================\n")
##
## ============================================================
cat("KRITERIA PEMILIHAN MODEL\n")
## KRITERIA PEMILIHAN MODEL
cat("============================================================\n")
## ============================================================
# AIC
aic_value <- AIC(model)
# BIC
bic_value <- BIC(model)
# AICc
n_obs <- nobs(model)
k <- length(coef(model))
aicc_value <- AIC(model) +
(2 * k * (k + 1)) / (n_obs - k - 1)
cat("\nAIC :", round(aic_value, 4), "\n")
##
## AIC : 1475.185
cat("BIC :", round(bic_value, 4), "\n")
## BIC : 1492.832
cat("AICc :", round(aicc_value, 4), "\n")
## AICc : 1475.347
# ------------------------------------------------------------
# 7. UJI ASUMSI RESIDUAL
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
# Ekstrak residual
residuals <- residuals(model)
fitted_values <- fitted(model)
standardized_resid <- rstandard(model)
# --- 7.1. UJI NORMALITAS RESIDUAL ---
cat("\n--- 6.1. UJI NORMALITAS RESIDUAL ---\n")
##
## --- 6.1. UJI NORMALITAS RESIDUAL ---
# Shapiro-Wilk Test
shapiro_test <- shapiro.test(residuals)
cat("\nShapiro-Wilk Test:\n")
##
## Shapiro-Wilk Test:
cat("W =", round(shapiro_test$statistic, 4), "\n")
## W = 0.9926
cat("p-value =", round(shapiro_test$p.value, 4), "\n")
## p-value = 0.2461
cat("Kesimpulan:", ifelse(shapiro_test$p.value > 0.05,
"Residual berdistribusi normal",
"Residual TIDAK berdistribusi normal"), "\n")
## Kesimpulan: Residual berdistribusi normal
# --- 7.2. UJI HETEROSKEDASTISITAS ---
cat("\n--- 6.2. UJI HETEROSKESASTISITAS ---\n")
##
## --- 6.2. UJI HETEROSKESASTISITAS ---
bp_test <- bptest(model)
glejser_test <- bptest(model, ~ fitted(model))
cat("Breusch-Pagan p-value:", round(bp_test$p.value, 4), "\n")
## Breusch-Pagan p-value: 0.1011
cat("Glejser p-value:", round(glejser_test$p.value, 4), "\n")
## Glejser p-value: 0.5647
# Kesimpulan
if (bp_test$p.value > 0.05 & glejser_test$p.value > 0.05) {
cat("Kesimpulan: Tidak ada heteroskedastisitas\n")
} else {
cat("Kesimpulan: Ada heteroskedastisitas\n")
}
## Kesimpulan: Tidak ada heteroskedastisitas
# --- 7.3. UJI AUTOKORELASI ---
cat("\n--- 6.3. UJI AUTOKORELASI ---\n")
##
## --- 6.3. UJI AUTOKORELASI ---
# Durbin-Watson Test
dw_test <- dwtest(model)
cat("\nDurbin-Watson Test:\n")
##
## Durbin-Watson Test:
cat("DW =", round(dw_test$statistic, 4), "\n")
## DW = 1.7945
cat("p-value =", round(dw_test$p.value, 4), "\n")
## p-value = 0.0469
cat("Kesimpulan:",
ifelse(dw_test$statistic < 1.5, "Ada autokorelasi positif",
ifelse(dw_test$statistic > 2.5, "Ada autokorelasi negatif",
"Tidak ada autokorelasi")), "\n")
## Kesimpulan: Tidak ada autokorelasi
# --- 7.4. UJI MULTIKOLINEARITAS ---
cat("\n--- 6.4. UJI MULTIKOLINEARITAS ---\n")
##
## --- 6.4. UJI MULTIKOLINEARITAS ---
# VIF (Variance Inflation Factor)
vif_values <- vif(model)
cat("\nVIF Values:\n")
##
## VIF Values:
print(round(vif_values, 4))
## X1 X2 X3
## 5.0377 10.3624 9.2847
cat("\nInterpretasi VIF:\n")
##
## Interpretasi VIF:
for (i in 1:length(vif_values)) {
vif_val <- vif_values[i]
var_name <- names(vif_values)[i]
status <- ifelse(vif_val < 5, "Tidak ada multikolinearitas",
ifelse(vif_val < 10, "Multikolinearitas moderat",
"Multikolinearitas serius"))
cat(sprintf("%-10s: VIF = %7.4f → %s\n", var_name, vif_val, status))
}
## X1 : VIF = 5.0377 → Multikolinearitas moderat
## X2 : VIF = 10.3624 → Multikolinearitas serius
## X3 : VIF = 9.2847 → Multikolinearitas moderat
# Tolerance
tolerance <- 1 / vif_values
cat("\nTolerance Values:\n")
##
## Tolerance Values:
print(round(tolerance, 4))
## X1 X2 X3
## 0.1985 0.0965 0.1077
# Condition Number
cat("\nCondition Number:", round(kappa(model), 4), "\n")
##
## Condition Number: 3859.019
cat("(Condition Number > 30 mengindikasikan multikolinearitas)\n")
## (Condition Number > 30 mengindikasikan multikolinearitas)
# --- 7.5. UJI LINEARITAS ---
cat("\n--- 6.5. UJI LINEARITAS ---\n")
##
## --- 6.5. UJI LINEARITAS ---
# Ramsey RESET Test
reset_test <- resettest(model, power = 2:3, type = "fitted")
cat("\nRamsey RESET Test:\n")
##
## Ramsey RESET Test:
cat("F =", round(reset_test$statistic, 4), "\n")
## F = 4.723
cat("p-value =", round(reset_test$p.value, 4), "\n")
## p-value = 0.0097
cat("Kesimpulan:", ifelse(reset_test$p.value > 0.05,
"Model linear (spesifikasi benar)",
"Model TIDAK linear (spesifikasi salah)"), "\n")
## Kesimpulan: Model TIDAK linear (spesifikasi salah)
# ------------------------------------------------------------
# 8. DIAGNOSTIK OUTLIER DAN INFLUENTIAL POINTS
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
# Cook's Distance
cooks_d <- cooks.distance(model)
cat("\nCook's Distance:\n")
##
## Cook's Distance:
cat("Nilai maksimum:", round(max(cooks_d), 4), "\n")
## Nilai maksimum: 0.4933
cat("Jumlah observasi dengan Cook's D > 1:", sum(cooks_d > 1), "\n")
## Jumlah observasi dengan Cook's D > 1: 0
# Leverage (Hat Values)
leverage <- hatvalues(model)
cat("\nLeverage (Hat Values):\n")
##
## Leverage (Hat Values):
cat("Nilai maksimum:", round(max(leverage), 4), "\n")
## Nilai maksimum: 0.1905
cat("Threshold 2(k+1)/n:", round(2 * (length(coef(model))) / n, 4), "\n")
## Threshold 2(k+1)/n: 0.0317
cat("Jumlah observasi dengan leverage > threshold:",
sum(leverage > 2 * length(coef(model)) / n), "\n")
## Jumlah observasi dengan leverage > threshold: 14
# Studentized Residuals
std_resid <- rstudent(model)
cat("\nStudentized Residuals:\n")
##
## Studentized Residuals:
cat("Nilai maksimum absolut:", round(max(abs(std_resid)), 4), "\n")
## Nilai maksimum absolut: 2.9405
cat("Jumlah observasi dengan |std_resid| > 3:",
sum(abs(std_resid) > 3), "\n")
## Jumlah observasi dengan |std_resid| > 3: 0
# DFBETAS
dfbetas_val <- dfbetas(model)
cat("\nDFBETAS:\n")
##
## DFBETAS:
cat("Jumlah observasi dengan |DFBETAS| > 1:",
sum(abs(dfbetas_val) > 1), "\n")
## Jumlah observasi dengan |DFBETAS| > 1: 0
# ------------------------------------------------------------
# 9. VISUALISASI DIAGNOSTIK
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
# Set par untuk 2x2 plot
par(mfrow = c(2, 2))
# Plot 1: Residual vs Fitted
plot(fitted_values, residuals,
xlab = "Fitted Values", ylab = "Residuals",
main = "Residual vs Fitted",
pch = 19, col = "steelblue")
abline(h = 0, col = "red", lty = 2)
lines(lowess(fitted_values, residuals), col = "green", lwd = 2)
# Plot 2: QQ-Plot
qqnorm(residuals, pch = 19, col = "steelblue",
main = "Normal Q-Q Plot")
qqline(residuals, col = "red", lwd = 2)
# Plot 3: Scale-Location
plot(fitted_values, sqrt(abs(standardized_resid)),
xlab = "Fitted Values", ylab = "√|Standardized Residuals|",
main = "Scale-Location",
pch = 19, col = "steelblue")
lines(lowess(fitted_values, sqrt(abs(standardized_resid))),
col = "green", lwd = 2)
# Plot 4: Residual vs Leverage
plot(leverage, standardized_resid,
xlab = "Leverage", ylab = "Standardized Residuals",
main = "Residual vs Leverage",
pch = 19, col = "steelblue")
abline(h = 0, col = "red", lty = 2)
abline(v = 2 * length(coef(model)) / n, col = "red", lty = 2)
# Reset par
par(mfrow = c(1, 1))
# ------------------------------------------------------------
# 10. UJI BOX-COX UNTUK TRANSFORMASI
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat("9. UJI BOX-COX\n")
## 9. UJI BOX-COX
cat(rep("=", 60) %+% "\n")
# Transformasi agar Y bernilai positif
Y_bc <- Y + 1
# Model dengan Y yang sudah ditransformasi
model_bc <- lm(Y_bc ~ X1 + X2 + X3)
bc <- boxcox(model_bc, lambda = seq(-2, 2, 0.1))
lambda_opt <- bc$x[which.max(bc$y)]
cat("\nOptimal Lambda:", round(lambda_opt, 4), "\n")
##
## Optimal Lambda: 1.0303
cat("Interpretasi:\n")
## Interpretasi:
cat("- Jika λ ≈ 1: Tidak perlu transformasi\n")
## - Jika λ ≈ 1: Tidak perlu transformasi
cat("- Jika λ ≈ 0: Gunakan log transformation\n")
## - Jika λ ≈ 0: Gunakan log transformation
cat("- Jika λ ≈ 0.5: Gunakan square root transformation\n")
## - Jika λ ≈ 0.5: Gunakan square root transformation
# ------------------------------------------------------------
# 11. PEMILIHAN MODEL (STEPWISE)
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(" PEMILIHAN MODEL (STEPWISE)\n")
## PEMILIHAN MODEL (STEPWISE)
cat(rep("=", 60) %+% "\n")
# Stepwise selection
step_model <- step(model, direction = "both", trace = 1)
## Start: AIC=758.04
## Y ~ X1 + X2 + X3
##
## Df Sum of Sq RSS AIC
## - X3 1 0 4943.2 756.04
## <none> 4943.2 758.04
## - X2 1 456 5399.2 778.27
## - X1 1 5675 10618.3 948.71
##
## Step: AIC=756.04
## Y ~ X1 + X2
##
## Df Sum of Sq RSS AIC
## <none> 4943.2 756.04
## + X3 1 0.0 4943.2 758.04
## - X2 1 1004.2 5947.5 800.65
## - X1 1 6042.7 10986.0 955.29
cat("\nModel Terbaik (Stepwise):\n")
##
## Model Terbaik (Stepwise):
print(summary(step_model))
##
## Call:
## lm(formula = Y ~ X1 + X2, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11.6459 -3.2071 -0.0299 3.1664 10.5064
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -45.95237 2.60501 -17.640 < 2e-16 ***
## X1 0.98950 0.05672 17.447 < 2e-16 ***
## X2 -0.14800 0.02081 -7.112 1.21e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.456 on 249 degrees of freedom
## Multiple R-squared: 0.7188, Adjusted R-squared: 0.7165
## F-statistic: 318.2 on 2 and 249 DF, p-value: < 2.2e-16
# ------------------------------------------------------------
# STEPWISE BACKWARD
# ------------------------------------------------------------
cat("\n============================================================\n")
##
## ============================================================
cat("STEPWISE BACKWARD\n")
## STEPWISE BACKWARD
cat("============================================================\n")
## ============================================================
model_backward <- step(
model,
direction = "backward",
trace = 1
)
## Start: AIC=758.04
## Y ~ X1 + X2 + X3
##
## Df Sum of Sq RSS AIC
## - X3 1 0 4943.2 756.04
## <none> 4943.2 758.04
## - X2 1 456 5399.2 778.27
## - X1 1 5675 10618.3 948.71
##
## Step: AIC=756.04
## Y ~ X1 + X2
##
## Df Sum of Sq RSS AIC
## <none> 4943.2 756.04
## - X2 1 1004.2 5947.5 800.65
## - X1 1 6042.7 10986.0 955.29
cat("\nModel Akhir Backward:\n")
##
## Model Akhir Backward:
print(formula(model_backward))
## Y ~ X1 + X2
cat("\nRingkasan Model Backward:\n")
##
## Ringkasan Model Backward:
print(summary(model_backward))
##
## Call:
## lm(formula = Y ~ X1 + X2, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11.6459 -3.2071 -0.0299 3.1664 10.5064
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -45.95237 2.60501 -17.640 < 2e-16 ***
## X1 0.98950 0.05672 17.447 < 2e-16 ***
## X2 -0.14800 0.02081 -7.112 1.21e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.456 on 249 degrees of freedom
## Multiple R-squared: 0.7188, Adjusted R-squared: 0.7165
## F-statistic: 318.2 on 2 and 249 DF, p-value: < 2.2e-16
# ------------------------------------------------------------
# STEPWISE FORWARD
# ------------------------------------------------------------
cat("\n============================================================\n")
##
## ============================================================
cat("STEPWISE FORWARD\n")
## STEPWISE FORWARD
cat("============================================================\n")
## ============================================================
model_null <- lm(Y ~ 1, data = data)
model_forward <- step(
model_null,
scope = formula(model),
direction = "forward",
trace = 1
)
## Start: AIC=1071.75
## Y ~ 1
##
## Df Sum of Sq RSS AIC
## + X1 1 11631.5 5947.5 800.65
## + X3 1 6871.2 10707.8 948.82
## + X2 1 6593.0 10986.0 955.29
## <none> 17579.0 1071.75
##
## Step: AIC=800.65
## Y ~ X1
##
## Df Sum of Sq RSS AIC
## + X2 1 1004.22 4943.2 756.04
## + X3 1 548.24 5399.2 778.27
## <none> 5947.5 800.65
##
## Step: AIC=756.04
## Y ~ X1 + X2
##
## Df Sum of Sq RSS AIC
## <none> 4943.2 756.04
## + X3 1 0.0052896 4943.2 758.04
cat("\nModel Akhir Forward:\n")
##
## Model Akhir Forward:
print(formula(model_forward))
## Y ~ X1 + X2
cat("\nRingkasan Model Forward:\n")
##
## Ringkasan Model Forward:
print(summary(model_forward))
##
## Call:
## lm(formula = Y ~ X1 + X2, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11.6459 -3.2071 -0.0299 3.1664 10.5064
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -45.95237 2.60501 -17.640 < 2e-16 ***
## X1 0.98950 0.05672 17.447 < 2e-16 ***
## X2 -0.14800 0.02081 -7.112 1.21e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.456 on 249 degrees of freedom
## Multiple R-squared: 0.7188, Adjusted R-squared: 0.7165
## F-statistic: 318.2 on 2 and 249 DF, p-value: < 2.2e-16
# ------------------------------------------------------------
# 12. RINGKASAN HASIL
# ------------------------------------------------------------
cat("\n" %+% rep("=", 60) %+% "\n")
cat(rep("=", 60) %+% "\n")
cat("\n--- MODEL AKHIR ---\n")
##
## --- MODEL AKHIR ---
cat("Persamaan Regresi:\n")
## Persamaan Regresi:
cat(sprintf("Y = %.4f + %.4f*X1 + %.4f*X2 + %.4f*X3 \n",
coef(model)[1], coef(model)[2], coef(model)[3],
coef(model)[4]))
## Y = -45.8455 + 0.9897*X1 + -0.1476*X2 + -0.0020*X3
cat("\n--- UJI SIGNIFIKANSI ---\n")
##
## --- UJI SIGNIFIKANSI ---
cat("Uji F (Simultan): p-value =", format(f_pvalue, scientific = TRUE), "\n")
## Uji F (Simultan): p-value = 5.101077e-68
cat("Uji t (Parsial):\n")
## Uji t (Parsial):
for (i in 2:nrow(coef_df)) {
cat(sprintf(" %s: p-value = %.4f\n",
coef_df$Variabel[i], coef_df$p_value[i]))
}
## X1: p-value = 0.0000
## X2: p-value = 0.0000
## X3: p-value = 0.9870
cat("\n--- UJI ASUMSI RESIDUAL ---\n")
##
## --- UJI ASUMSI RESIDUAL ---
cat("Normalitas (Shapiro-Wilk): p-value =",
round(shapiro_test$p.value, 4), "\n")
## Normalitas (Shapiro-Wilk): p-value = 0.2461
cat("Heteroskedastisitas (Breusch-Pagan): p-value =",
round(bp_test$p.value, 4), "\n")
## Heteroskedastisitas (Breusch-Pagan): p-value = 0.1011
cat("Autokorelasi (Durbin-Watson): DW =",
round(dw_test$statistic, 4), "\n")
## Autokorelasi (Durbin-Watson): DW = 1.7945
cat("Multikolinearitas (VIF max):",
round(max(vif_values), 4), "\n")
## Multikolinearitas (VIF max): 10.3624
cat("Linearitas (Ramsey RESET): p-value =",
round(reset_test$p.value, 4), "\n")
## Linearitas (Ramsey RESET): p-value = 0.0097
cat("\n--- KESIMPULAN ---\n")
##
## --- KESIMPULAN ---
cat("R-squared:", round(r_squared, 4), "\n")
## R-squared: 0.7188
cat("Adjusted R-squared:", round(adj_r_squared, 4), "\n")
## Adjusted R-squared: 0.7154
cat("\n" %+% rep("=", 60) %+% "\n")
cat("ANALISIS SELESAI\n")
## ANALISIS SELESAI
cat(rep("=", 60) %+% "\n")