# 1. MEMUAT PACKAGE DAN DATA
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.5.2
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.5.2
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.5.2
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(car)
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.5.2
# Membaca data (Palmer Penguins)
penguins <- read.csv("C:/Analisis Regresi/penguins_size.csv")

# Menghapus missing value
data <- na.omit(penguins)

# Memilih variabel
# X : flipper_length_mm
# Y : body_mass_g
# 2. ANALISIS DESKRIPTIF DAN VISUALISASI
print("Statistik Deskriptif:")
## [1] "Statistik Deskriptif:"
summary(data[, c("flipper_length_mm", "body_mass_g")])
##  flipper_length_mm  body_mass_g  
##  Min.   :172       Min.   :2700  
##  1st Qu.:190       1st Qu.:3550  
##  Median :197       Median :4050  
##  Mean   :201       Mean   :4209  
##  3rd Qu.:213       3rd Qu.:4794  
##  Max.   :231       Max.   :6300
# Scatter plot
ggplot(data, aes(x = flipper_length_mm, y = body_mass_g)) +
  geom_point(color = "blue", size = 3) +
  labs(title = "Hubungan Panjang Flipper dan Berat Badan Penguin",
       x = "Panjang Flipper (mm)",
       y = "Berat Badan (gram)") +
  theme_minimal()

# Korelasi
cor_test <- cor.test(data$flipper_length_mm, data$body_mass_g)
print(paste("Korelasi Pearson:", round(cor_test$estimate, 4)))
## [1] "Korelasi Pearson: 0.8732"
print(paste("p-value korelasi:", round(cor_test$p.value, 4)))
## [1] "p-value korelasi: 0"
# 3. MEMBANGUN MODEL REGRESI
model <- lm(body_mass_g ~ flipper_length_mm, data = data)
print("Ringkasan Model Regresi:")
## [1] "Ringkasan Model Regresi:"
summary(model)
## 
## Call:
## lm(formula = body_mass_g ~ flipper_length_mm, data = data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1057.23  -259.15   -11.13   242.86  1293.73 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       -5865.818    309.340  -18.96   <2e-16 ***
## flipper_length_mm    50.120      1.535   32.65   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 392.8 on 332 degrees of freedom
## Multiple R-squared:  0.7625, Adjusted R-squared:  0.7618 
## F-statistic:  1066 on 1 and 332 DF,  p-value: < 2.2e-16
# 4. UJI ASUMSI REGRESI LINEAR
cat("\n=== UJI ASUMSI REGRESI LINEAR ===\n")
## 
## === UJI ASUMSI REGRESI LINEAR ===
# 4.1 Normalitas Residual
shapiro_test <- shapiro.test(residuals(model))
cat("1. UJI NORMALITAS (Shapiro-Wilk):\n")
## 1. UJI NORMALITAS (Shapiro-Wilk):
cat("   Statistik W =", round(shapiro_test$statistic, 4), "\n")
##    Statistik W = 0.9922
cat("   p-value =", round(shapiro_test$p.value, 4), "\n")
##    p-value = 0.0792
if(shapiro_test$p.value > 0.05) {
  cat("   Keputusan: Residual berdistribusi normal\n")
} else {
  cat("   Keputusan: Residual tidak normal\n")
}
##    Keputusan: Residual berdistribusi normal
# Q-Q Plot
qqnorm(residuals(model), main = "Q-Q Plot Residual")
qqline(residuals(model), col = "red")

# 4.2 Homoskedastisitas
bp_test <- bptest(model)
cat("\n2. UJI HOMOSKEDASTISITAS (Breusch-Pagan):\n")
## 
## 2. UJI HOMOSKEDASTISITAS (Breusch-Pagan):
cat("   Statistik LM =", round(bp_test$statistic, 4), "\n")
##    Statistik LM = 2.897
cat("   p-value =", round(bp_test$p.value, 4), "\n")
##    p-value = 0.0887
if(bp_test$p.value > 0.05) {
  cat("   Keputusan: Varian residual homogen\n")
} else {
  cat("   Keputusan: Ada heteroskedastisitas\n")
}
##    Keputusan: Varian residual homogen
# Plot Residual vs Fitted
plot(fitted(model), residuals(model),
     main = "Residual vs Fitted Values",
     xlab = "Fitted Values", ylab = "Residuals",
     pch = 19, col = "blue")
abline(h = 0, col = "red", lty = 2)

# 4.3 Tidak ada Autokorelasi
dw_test <- dwtest(model)
cat("\n3. UJI AUTOKORELASI (Durbin-Watson):\n")
## 
## 3. UJI AUTOKORELASI (Durbin-Watson):
cat("   Statistik DW =", round(dw_test$statistic, 4), "\n")
##    Statistik DW = 2.1125
cat("   p-value =", round(dw_test$p.value, 4), "\n")
##    p-value = 0.8392
if(dw_test$p.value > 0.05) {
  cat("   Keputusan: Tidak ada autokorelasi\n")
} else {
  cat("   Keputusan: Ada autokorelasi\n")
}
##    Keputusan: Tidak ada autokorelasi
# 5. INTERPRETASI KOEFISIEN
cat("\n=== INTERPRETASI KOEFISIEN ===\n")
## 
## === INTERPRETASI KOEFISIEN ===
intercept <- coef(model)[1]
slope <- coef(model)[2]

cat("Persamaan Regresi: Body_Mass =", 
    round(intercept, 2), "+", 
    round(slope, 2), "* Flipper_Length\n")
## Persamaan Regresi: Body_Mass = -5865.82 + 50.12 * Flipper_Length
cat("\nInterpretasi:\n")
## 
## Interpretasi:
cat("1. Intercept (β0): Berat badan saat panjang flipper = 0 mm\n")
## 1. Intercept (β0): Berat badan saat panjang flipper = 0 mm
cat("2. Slope (β1): Setiap kenaikan 1 mm flipper, berat badan meningkat",
    round(slope, 2), "gram\n")
## 2. Slope (β1): Setiap kenaikan 1 mm flipper, berat badan meningkat 50.12 gram
# 6. ESTIMASI PARAMETER DAN INFERENSI
cat("\n=== ESTIMASI PARAMETER ===\n")
## 
## === ESTIMASI PARAMETER ===
conf_int <- confint(model, level = 0.95)
cat("Interval Kepercayaan 95%:\n")
## Interval Kepercayaan 95%:
cat("   Intercept: [", round(conf_int[1,1], 2), ", ",
    round(conf_int[1,2], 2), "]\n")
##    Intercept: [ -6474.33 ,  -5257.3 ]
cat("   Slope:     [", round(conf_int[2,1], 2), ", ",
    round(conf_int[2,2], 2), "]\n")
##    Slope:     [ 47.1 ,  53.14 ]
# Uji hipotesis slope
summary_model <- summary(model)
slope_pvalue <- summary_model$coefficients[2,4]
cat("\nUji Hipotesis:\n")
## 
## Uji Hipotesis:
cat("H0: β1 = 0 (tidak ada pengaruh)\n")
## H0: β1 = 0 (tidak ada pengaruh)
cat("H1: β1 ≠ 0 (ada pengaruh)\n")
## H1: β1 ≠ 0 (ada pengaruh)
cat("p-value =", round(slope_pvalue, 6), "\n")
## p-value = 0
if(slope_pvalue < 0.05) {
  cat("Keputusan: Tolak H0\n")
} else {
  cat("Keputusan: Gagal tolak H0\n")
}
## Keputusan: Tolak H0
# 7. KOEFISIEN DETERMINASI
r_squared <- summary_model$r.squared
cat("\nKoefisien Determinasi (R²):", round(r_squared, 4), "\n")
## 
## Koefisien Determinasi (R²): 0.7625
cat("Artinya", round(r_squared*100, 2),
    "% variasi berat badan dijelaskan oleh panjang flipper\n")
## Artinya 76.25 % variasi berat badan dijelaskan oleh panjang flipper
# 8. VISUALISASI MODEL
ggplot(data, aes(x = flipper_length_mm, y = body_mass_g)) +
  geom_point(color = "blue", size = 3) +
  geom_smooth(method = "lm", se = TRUE,
              color = "red", fill = "pink") +
  labs(title = "Garis Regresi Linear",
       x = "Panjang Flipper (mm)",
       y = "Berat Badan (gram)") +
  theme_minimal()
## `geom_smooth()` using formula = 'y ~ x'

# 9. RINGKASAN ANALISIS
cat("\n=== RINGKASAN ANALISIS ===\n")
## 
## === RINGKASAN ANALISIS ===
cat("1. Model regresi signifikan\n")
## 1. Model regresi signifikan
cat("2. Hubungan linear positif\n")
## 2. Hubungan linear positif
cat("3. Asumsi regresi terpenuhi\n")
## 3. Asumsi regresi terpenuhi
cat("4. Panjang flipper berpengaruh terhadap berat badan\n")
## 4. Panjang flipper berpengaruh terhadap berat badan