library(splines)

df_airquality <- datasets::airquality

df_airquality$Ozone[is.na(df_airquality$Ozone)] <- median(df_airquality$Ozone, na.rm = TRUE)
df_airquality$Solar.R[is.na(df_airquality$Solar.R)] <- median(df_airquality$Solar.R, na.rm = TRUE)
MSE <- function(pred, actual) {
  mean((pred - actual)^2, na.rm = TRUE)
}

1. Dengan menggunakan model dan syntax yang sama, lakukan perbandingan antara variabel Ozone dengan Temp.

mod_linear_temp <- lm(Ozone ~ Temp, data = df_airquality)
summary(mod_linear_temp)
## 
## Call:
## lm(formula = Ozone ~ Temp, data = df_airquality)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -36.464 -16.399  -3.866  10.314 122.691 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -104.0802    15.6655  -6.644 5.21e-10 ***
## Temp           1.8443     0.1997   9.236  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 23.3 on 151 degrees of freedom
## Multiple R-squared:  0.361,  Adjusted R-squared:  0.3568 
## F-statistic: 85.31 on 1 and 151 DF,  p-value: < 2.2e-16
pred_linear_temp <- predict(mod_linear_temp)

# Menghitung performa
MSE_linear_temp <- mean((pred_linear_temp - df_airquality$Ozone)^2, na.rm = TRUE)
AIC_linear_temp <- AIC(mod_linear_temp)
AdjR2_linear_temp <- summary(mod_linear_temp)$adj.r.squared

data.frame(
  Model = "Linear",
  MSE = MSE_linear_temp,
  AIC = AIC_linear_temp,
  Adj_R2 = AdjR2_linear_temp
)
##    Model      MSE      AIC    Adj_R2
## 1 Linear 535.8932 1401.637 0.3567682

Hasil analisis regresi linier menunjukkan bahwa variabel suhu (Temp) berpengaruh positif dan signifikan terhadap kadar Ozone. Koefisien regresi sebesar 1.8443 berarti setiap kenaikan 1 derajat suhu akan meningkatkan kadar Ozone sekitar 1.84 unit. Nilai p < 0.001 menegaskan hubungan ini sangat signifikan. Nilai Adjusted R² sebesar 0.357 mengindikasikan bahwa sekitar 35,7% variasi Ozone dapat dijelaskan oleh Temp, sementara sisanya dipengaruhi faktor lain. MSE sebesar 535.89 masih cukup besar, sehingga prediksi belum sepenuhnya akurat, tetapi model linear ini sudah cukup memberikan gambaran awal bahwa suhu memainkan peran penting dalam pembentukan Ozone.

2. Ulangi analisis untuk hubungan antara Ozone dengan Wind.

# Membuat model fungsi tangga dengan membagi Wind jadi 5 kategori
mod_tangga_wind <- lm(Ozone ~ cut(Wind, 5), data = df_airquality)
summary(mod_tangga_wind)
## 
## Call:
## lm(formula = Ozone ~ cut(Wind, 5), data = df_airquality)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -55.731 -14.807  -4.952   5.696  80.769 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               87.231      6.588  13.241  < 2e-16 ***
## cut(Wind, 5)(5.5,9.3]    -40.924      7.300  -5.606 9.86e-08 ***
## cut(Wind, 5)(9.3,13.1]   -58.279      7.365  -7.913 5.43e-13 ***
## cut(Wind, 5)(13.1,16.9]  -61.427      7.972  -7.706 1.74e-12 ***
## cut(Wind, 5)(16.9,20.7]  -70.231     15.214  -4.616 8.41e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 23.75 on 148 degrees of freedom
## Multiple R-squared:  0.3493, Adjusted R-squared:  0.3317 
## F-statistic: 19.86 on 4 and 148 DF,  p-value: 4.178e-13
# Memprediksi dari model
pred_tangga_wind <- predict(mod_tangga_wind)

# Menghitung performa model
MSE_tangga_wind <- mean((pred_tangga_wind - df_airquality$Ozone)^2, na.rm = TRUE)
AIC_tangga_wind <- AIC(mod_tangga_wind)
AdjR2_tangga_wind <- summary(mod_tangga_wind)$adj.r.squared

# Membuat ringkasan data frame
data.frame(
  Model = "Tangga (Wind)",
  MSE = MSE_tangga_wind,
  AIC = AIC_tangga_wind,
  Adj_R2 = AdjR2_tangga_wind
)
##           Model      MSE      AIC   Adj_R2
## 1 Tangga (Wind) 545.7385 1410.423 0.331673

Model fungsi tangga yang membagi kecepatan angin (Wind) menjadi lima kategori menunjukkan bahwa rata-rata kadar Ozone menurun signifikan seiring bertambahnya kategori kecepatan angin. Koefisien negatif pada semua interval Wind menegaskan bahwa dibandingkan kelompok Wind terendah, kadar Ozone cenderung lebih rendah pada kelompok angin sedang hingga tinggi. Nilai Adjusted R² sebesar 0.332 menunjukkan model ini mampu menjelaskan sekitar 33% variasi Ozone, sedikit lebih rendah dibanding model linear sederhana. MSE sebesar 545.74 juga relatif tinggi, sehingga meskipun hubungan Wind dengan Ozone jelas signifikan, pendekatan fungsi tangga kurang efisien dibanding model lain karena menghasilkan prediksi “meloncat” antar kategori, bukan pola yang halus.

3. Bandingkan hasil dari masing-masing model berdasarkan nilai MSE, AIC, dan Adjusted R-squared.

mod_linear_wind <- lm(Ozone ~ Wind, data = df_airquality)
pred_linear_wind <- predict(mod_linear_wind)

MSE_linear_wind <- mean((pred_linear_wind - df_airquality$Ozone)^2, na.rm = TRUE)
AIC_linear_wind <- AIC(mod_linear_wind)
AdjR2_linear_wind <- summary(mod_linear_wind)$adj.r.squared

# Tangga Ozone ~ Wind
mod_tangga_wind <- lm(Ozone ~ cut(Wind, 5), data = df_airquality)
pred_tangga_wind <- predict(mod_tangga_wind)

MSE_tangga_wind <- mean((pred_tangga_wind - df_airquality$Ozone)^2, na.rm = TRUE)
AIC_tangga_wind <- AIC(mod_tangga_wind)
AdjR2_tangga_wind <- summary(mod_tangga_wind)$adj.r.squared

# Cubic Spline Ozone ~ Wind

knots_wind <- quantile(df_airquality$Wind, probs = c(0.25, 0.5, 0.75))
mod_spline_wind <- lm(Ozone ~ bs(Wind, knots = as.numeric(knots_wind)), data = df_airquality)
pred_spline_wind <- predict(mod_spline_wind)

MSE_spline_wind <- mean((pred_spline_wind - df_airquality$Ozone)^2, na.rm = TRUE)
AIC_spline_wind <- AIC(mod_spline_wind)
AdjR2_spline_wind <- summary(mod_spline_wind)$adj.r.squared

# Natural Spline Ozone ~ Wind

mod_nspline_wind <- lm(Ozone ~ ns(Wind, knots = as.numeric(knots_wind)), data = df_airquality)
pred_nspline_wind <- predict(mod_nspline_wind)

MSE_nspline_wind <- mean((pred_nspline_wind - df_airquality$Ozone)^2, na.rm = TRUE)
AIC_nspline_wind <- AIC(mod_nspline_wind)
AdjR2_nspline_wind <- summary(mod_nspline_wind)$adj.r.squared

# Tabel Perbandingan Semua Model Wind

compare_wind <- data.frame(
  Model = c("Linear", "Tangga", "Spline", "Natural Spline"),
  MSE   = c(MSE_linear_wind, MSE_tangga_wind, MSE_spline_wind, MSE_nspline_wind),
  AIC   = c(AIC_linear_wind, AIC_tangga_wind, AIC_spline_wind, AIC_nspline_wind),
  Adj_R2= c(AdjR2_linear_wind, AdjR2_tangga_wind, AdjR2_spline_wind, AdjR2_nspline_wind)
)

print("Perbandingan Ozone ~ Wind")
## [1] "Perbandingan Ozone ~ Wind"
print(compare_wind)
##            Model      MSE      AIC    Adj_R2
## 1         Linear 601.3750 1419.276 0.2781705
## 2         Tangga 545.7385 1410.423 0.3316730
## 3         Spline 480.4891 1394.940 0.4035187
## 4 Natural Spline 513.2467 1401.031 0.3714633

Hasil perbandingan menunjukkan bahwa model Spline memiliki kinerja terbaik dengan MSE paling rendah (480.49), AIC terendah (1394.94), dan Adjusted R² tertinggi (0.404). Hal ini berarti spline mampu menangkap pola non-linear antara Ozone dan Wind dengan lebih baik dibandingkan model lainnya. Model Natural Spline juga cukup baik dengan nilai MSE dan Adjusted R² yang kompetitif, serta lebih stabil di ekor data. Sementara itu, model Tangga menghasilkan performa lebih baik daripada linear, namun prediksinya cenderung “loncat” antar kategori sehingga kurang halus. Model Linear adalah yang paling sederhana, namun memiliki MSE terbesar dan Adjusted R² paling rendah. Dengan demikian, Spline direkomendasikan untuk hubungan Ozone–Wind karena memberikan keseimbangan terbaik antara akurasi dan kemampuan menjelaskan variansi.

4. Berikan insight dan pendapatmu.

Model Natural Spline memberikan hasil yang cukup baik dalam menjelaskan hubungan Ozone dengan Wind. Dibandingkan model linear, natural spline berhasil menurunkan MSE (dari 601.37 menjadi 513.25) dan meningkatkan Adjusted R² (dari 0.278 menjadi 0.371). Artinya, natural spline mampu menangkap pola non-linear yang ada pada data, sehingga prediksi Ozone terhadap variasi kecepatan angin lebih akurat dan variansi yang dijelaskan lebih besar. Nilai AIC (1401.03) juga lebih rendah daripada linear, menandakan model lebih efisien meski sedikit lebih kompleks.

Menurut saya, natural spline ini cocok dipakai karena bisa menggambarkan perubahan kadar Ozone secara lebih halus tanpa menghasilkan pola “loncat-loncat” seperti pada fungsi tangga. Meski performanya masih sedikit di bawah cubic spline murni (yang punya Adj R² lebih tinggi), natural spline lebih stabil pada bagian ekor distribusi data, sehingga lebih aman dari risiko overfitting. Dengan demikian, untuk aplikasi nyata, natural spline bisa menjadi pilihan kompromi antara akurasi dan stabilitas model dalam memprediksi hubungan Ozone dengan Wind.