Tugas Mandiri PSD P2
Syifa Khairunnisa
2023-09-28
Library
library(glmnet)## Warning: package 'glmnet' was built under R version 4.2.3## Loading required package: Matrix## Loaded glmnet 4.1-8Data
Sumber data: https://www.kaggle.com/datasets/anuragbantu/new-delhi-air-quality
Peubah respon:
y = AQI (Air Quality Index)
Peubah penjelas:
x1 = CO
x2 = no2
x3 = o3
x4 = pm10
x5 = pm25
x6 = so2
library(rio)## Warning: package 'rio' was built under R version 4.2.3data <- import("https://raw.githubusercontent.com/syifaskn/PSD/main/Tugas%20Individu/DataAQI_Syifa.csv")
data## V1 AQI CO datetime no2 o3 pm10 pm25
## 1 0 30.2 198.6027 2022-10-21:18 0.04685717 55.78995 10.486722 5.637410
## 2 1 28.2 197.6013 2022-10-21:19 0.04645554 54.93164 10.719325 4.618169
## 3 2 26.6 198.6027 2022-10-21:20 0.04685717 54.64554 11.155578 3.520902
## 4 3 25.0 201.9405 2022-10-21:21 0.04819594 55.07469 11.116206 2.225919
## 5 4 26.0 205.2784 2022-10-21:22 0.04886533 55.78995 10.405250 1.979471
## 6 5 26.0 208.6163 2022-10-21:23 0.05020411 56.50520 9.402518 1.792677
## 7 6 26.0 208.6163 2022-10-22:00 0.05020411 55.78995 8.482822 1.752237
## 8 7 24.0 205.2784 2022-10-22:01 0.05020411 52.92892 7.586667 1.858459
## 9 8 23.0 200.2716 2022-10-22:02 0.04685717 50.06790 7.636587 2.208076
## 10 9 24.0 196.9338 2022-10-22:03 0.04417961 51.49841 8.487648 2.476681
## 11 10 25.0 198.6027 2022-10-22:04 0.04284084 53.64418 9.357371 2.585729
## 12 11 26.0 200.2716 2022-10-22:05 0.04049798 55.78995 10.338728 2.390081
## 13 12 27.0 203.6095 2022-10-22:06 0.04016329 59.36623 10.523314 2.086022
## 14 13 29.0 205.2784 2022-10-22:07 0.03882451 62.22725 10.635970 1.937047
## 15 14 29.0 205.2784 2022-10-22:08 0.03882451 62.94250 9.964803 1.835311
## 16 15 29.0 205.2784 2022-10-22:09 0.03882451 62.94250 9.214411 1.725189
## 17 16 29.0 203.6095 2022-10-22:10 0.03748573 62.22725 8.772117 1.640662
## 18 17 28.0 201.9405 2022-10-22:11 0.04250614 60.79674 8.753794 1.633373
## 19 18 27.0 200.2716 2022-10-22:12 0.05020411 59.36623 8.807467 1.636810
## 20 19 27.0 201.9405 2022-10-22:13 0.05622860 58.65097 8.774860 1.666231
## 21 20 27.0 201.9405 2022-10-22:14 0.05622860 58.65097 8.675813 1.690578
## 22 21 27.0 201.9405 2022-10-22:15 0.05488982 57.93572 8.915920 1.705802
## 23 22 27.0 200.2716 2022-10-22:16 0.05288166 57.93572 9.227945 1.732754
## 24 23 27.0 200.2716 2022-10-22:17 0.05221227 58.65097 9.491790 1.765387
## 25 24 27.0 200.2716 2022-10-22:18 0.05221227 59.36623 9.907828 1.827327
## 26 25 28.0 201.9405 2022-10-22:19 0.05221227 61.51199 10.095885 1.888416
## 27 26 29.0 203.6095 2022-10-22:20 0.05288166 63.65776 9.879244 1.892860
## 28 27 31.0 205.2784 2022-10-22:21 0.05221227 66.51878 9.503339 1.859564
## 29 28 31.0 206.9473 2022-10-22:22 0.05221227 67.94930 9.008138 1.829196
## 30 29 32.0 206.9473 2022-10-22:23 0.05154289 68.66455 8.729765 1.787945
## 31 30 32.0 205.2784 2022-10-23:00 0.05087350 68.66455 8.912433 1.796048
## 32 31 31.0 206.9473 2022-10-23:01 0.05488982 67.23404 9.337687 1.810690
## 33 32 31.0 206.9473 2022-10-23:02 0.06091431 66.51878 9.351754 1.826056
## 34 33 30.0 208.6163 2022-10-23:03 0.05555921 65.80353 8.640116 1.822618
## 35 34 30.0 206.9473 2022-10-23:04 0.04618778 64.37302 7.750329 1.819336
## 36 35 29.0 205.2784 2022-10-23:05 0.03815512 62.22725 6.815522 1.831419
## 37 36 27.0 203.6095 2022-10-23:06 0.03079185 59.36623 6.965347 1.929624
## 38 37 26.0 200.2716 2022-10-23:07 0.02577144 56.50520 7.022172 2.049507
## 39 38 25.0 200.2716 2022-10-23:08 0.02342858 55.07469 6.946332 2.152871
## 40 39 25.0 198.6027 2022-10-23:09 0.02376328 53.64418 6.834290 2.219349
## 41 40 24.0 196.9338 2022-10-23:10 0.02677552 52.92892 6.686746 2.248050
## 42 41 24.0 195.2648 2022-10-23:11 0.03146124 52.92892 7.593563 2.276788
## 43 42 25.0 195.2648 2022-10-23:12 0.03815512 53.64418 8.999838 2.334924
## 44 43 25.0 195.2648 2022-10-23:13 0.04351023 53.64418 11.022134 2.365020
## 45 44 25.0 196.9338 2022-10-23:14 0.04351023 54.35944 12.198974 2.366476
## 46 45 26.0 198.6027 2022-10-23:15 0.04484900 55.78995 11.969586 2.369681
## 47 46 26.0 198.6027 2022-10-23:16 0.04618778 57.22046 9.593794 2.339545
## 48 47 27.0 198.6027 2022-10-23:17 0.04551839 57.93572 7.763611 2.331504
## 49 48 27.0 198.6027 2022-10-23:18 0.04417961 58.65097 7.517936 2.324217
## 50 49 27.0 198.6027 2022-10-23:19 0.04484900 58.65097 8.041439 2.320275
## 51 50 27.0 200.2716 2022-10-23:20 0.04618778 59.36623 8.654108 2.283845
## 52 51 28.0 200.2716 2022-10-23:21 0.04752655 60.08148 8.514708 2.195546
## 53 52 28.0 200.2716 2022-10-23:22 0.04685717 60.08148 7.923108 2.066149
## 54 53 27.0 200.2716 2022-10-23:23 0.04618778 59.36623 7.097172 2.059171
## 55 54 27.0 198.6027 2022-10-24:00 0.04551839 58.65097 6.815780 2.072363
## 56 55 26.0 198.6027 2022-10-24:01 0.04752655 57.22046 7.661589 2.102422
## 57 56 25.0 198.6027 2022-10-24:02 0.04551839 54.35944 8.718614 2.230021
## 58 57 23.0 198.6027 2022-10-24:03 0.03715104 50.78316 8.457276 2.501037
## 59 58 22.0 196.9338 2022-10-24:04 0.02744491 47.92213 7.184991 2.750613
## 60 59 21.0 193.5959 2022-10-24:05 0.02125307 45.41874 6.950173 2.984422
## 61 60 20.0 193.5959 2022-10-24:06 0.01757144 43.98823 7.205668 3.065950
## 62 61 19.0 191.9270 2022-10-24:07 0.01573062 41.84246 7.898480 3.118563
## 63 62 19.0 191.9270 2022-10-24:08 0.01556327 41.48483 8.744758 3.090378
## 64 63 20.0 191.9270 2022-10-24:09 0.01807348 43.27297 9.623435 2.976335
## 65 64 21.0 191.9270 2022-10-24:10 0.02208981 45.41874 10.253539 2.846886
## 66 65 21.0 191.9270 2022-10-24:11 0.02577144 45.41874 10.703334 2.771884
## 67 66 21.0 191.9270 2022-10-24:12 0.03045716 45.41874 10.663057 2.730851
## 68 67 22.0 193.5959 2022-10-24:13 0.03547757 47.20688 10.423121 2.720472
## 69 68 24.0 195.2648 2022-10-24:14 0.03915921 52.21367 11.240661 2.713109
## 70 69 26.0 196.9338 2022-10-24:15 0.04417961 56.50520 12.098125 2.743044
## 71 70 27.0 198.6027 2022-10-24:16 0.04685717 59.36623 12.845977 2.780022
## 72 71 28.0 198.6027 2022-10-24:17 0.04752655 60.79674 12.583639 2.745875
## so2 timestamp_local timestamp_utc ts
## 1 0.3874302 2022-10-21 23:00:00 2022-10-21 18:00:00 1666375200
## 2 0.4097819 2022-10-22 00:00:00 2022-10-21 19:00:00 1666378800
## 3 0.4023313 2022-10-22 01:00:00 2022-10-21 20:00:00 1666382400
## 4 0.3762543 2022-10-22 02:00:00 2022-10-21 21:00:00 1666386000
## 5 0.3390014 2022-10-22 03:00:00 2022-10-21 22:00:00 1666389600
## 6 0.3129244 2022-10-22 04:00:00 2022-10-21 23:00:00 1666393200
## 7 0.3129244 2022-10-22 05:00:00 2022-10-22 00:00:00 1666396800
## 8 0.3427267 2022-10-22 06:00:00 2022-10-22 01:00:00 1666400400
## 9 0.3948808 2022-10-22 07:00:00 2022-10-22 02:00:00 1666404000
## 10 0.4284084 2022-10-22 08:00:00 2022-10-22 03:00:00 1666407600
## 11 0.4284084 2022-10-22 09:00:00 2022-10-22 04:00:00 1666411200
## 12 0.3986061 2022-10-22 10:00:00 2022-10-22 05:00:00 1666414800
## 13 0.3501773 2022-10-22 11:00:00 2022-10-22 06:00:00 1666418400
## 14 0.3166497 2022-10-22 12:00:00 2022-10-22 07:00:00 1666422000
## 15 0.2980232 2022-10-22 13:00:00 2022-10-22 08:00:00 1666425600
## 16 0.2868474 2022-10-22 14:00:00 2022-10-22 09:00:00 1666429200
## 17 0.2831221 2022-10-22 15:00:00 2022-10-22 10:00:00 1666432800
## 18 0.2942979 2022-10-22 16:00:00 2022-10-22 11:00:00 1666436400
## 19 0.3017485 2022-10-22 17:00:00 2022-10-22 12:00:00 1666440000
## 20 0.3054738 2022-10-22 18:00:00 2022-10-22 13:00:00 1666443600
## 21 0.3091991 2022-10-22 19:00:00 2022-10-22 14:00:00 1666447200
## 22 0.3129244 2022-10-22 20:00:00 2022-10-22 15:00:00 1666450800
## 23 0.3203750 2022-10-22 21:00:00 2022-10-22 16:00:00 1666454400
## 24 0.3241003 2022-10-22 22:00:00 2022-10-22 17:00:00 1666458000
## 25 0.3278256 2022-10-22 23:00:00 2022-10-22 18:00:00 1666461600
## 26 0.3278256 2022-10-23 00:00:00 2022-10-22 19:00:00 1666465200
## 27 0.3203750 2022-10-23 01:00:00 2022-10-22 20:00:00 1666468800
## 28 0.3166497 2022-10-23 02:00:00 2022-10-22 21:00:00 1666472400
## 29 0.3129244 2022-10-23 03:00:00 2022-10-22 22:00:00 1666476000
## 30 0.3166497 2022-10-23 04:00:00 2022-10-22 23:00:00 1666479600
## 31 0.3166497 2022-10-23 05:00:00 2022-10-23 00:00:00 1666483200
## 32 0.3166497 2022-10-23 06:00:00 2022-10-23 01:00:00 1666486800
## 33 0.3203750 2022-10-23 07:00:00 2022-10-23 02:00:00 1666490400
## 34 0.3203750 2022-10-23 08:00:00 2022-10-23 03:00:00 1666494000
## 35 0.3241003 2022-10-23 09:00:00 2022-10-23 04:00:00 1666497600
## 36 0.3278256 2022-10-23 10:00:00 2022-10-23 05:00:00 1666501200
## 37 0.3390014 2022-10-23 11:00:00 2022-10-23 06:00:00 1666504800
## 38 0.3501773 2022-10-23 12:00:00 2022-10-23 07:00:00 1666508400
## 39 0.3576279 2022-10-23 13:00:00 2022-10-23 08:00:00 1666512000
## 40 0.3650784 2022-10-23 14:00:00 2022-10-23 09:00:00 1666515600
## 41 0.3725290 2022-10-23 15:00:00 2022-10-23 10:00:00 1666519200
## 42 0.3762543 2022-10-23 16:00:00 2022-10-23 11:00:00 1666522800
## 43 0.3762543 2022-10-23 17:00:00 2022-10-23 12:00:00 1666526400
## 44 0.3725290 2022-10-23 18:00:00 2022-10-23 13:00:00 1666530000
## 45 0.3725290 2022-10-23 19:00:00 2022-10-23 14:00:00 1666533600
## 46 0.3725290 2022-10-23 20:00:00 2022-10-23 15:00:00 1666537200
## 47 0.3799796 2022-10-23 21:00:00 2022-10-23 16:00:00 1666540800
## 48 0.3911555 2022-10-23 22:00:00 2022-10-23 17:00:00 1666544400
## 49 0.3986061 2022-10-23 23:00:00 2022-10-23 18:00:00 1666548000
## 50 0.3986061 2022-10-24 00:00:00 2022-10-23 19:00:00 1666551600
## 51 0.3874302 2022-10-24 01:00:00 2022-10-23 20:00:00 1666555200
## 52 0.3762543 2022-10-24 02:00:00 2022-10-23 21:00:00 1666558800
## 53 0.3650784 2022-10-24 03:00:00 2022-10-23 22:00:00 1666562400
## 54 0.3688037 2022-10-24 04:00:00 2022-10-23 23:00:00 1666566000
## 55 0.3762543 2022-10-24 05:00:00 2022-10-24 00:00:00 1666569600
## 56 0.3837049 2022-10-24 06:00:00 2022-10-24 01:00:00 1666573200
## 57 0.4023313 2022-10-24 07:00:00 2022-10-24 02:00:00 1666576800
## 58 0.4284084 2022-10-24 08:00:00 2022-10-24 03:00:00 1666580400
## 59 0.4433095 2022-10-24 09:00:00 2022-10-24 04:00:00 1666584000
## 60 0.4544854 2022-10-24 10:00:00 2022-10-24 05:00:00 1666587600
## 61 0.4544854 2022-10-24 11:00:00 2022-10-24 06:00:00 1666591200
## 62 0.4507601 2022-10-24 12:00:00 2022-10-24 07:00:00 1666594800
## 63 0.4395843 2022-10-24 13:00:00 2022-10-24 08:00:00 1666598400
## 64 0.4209578 2022-10-24 14:00:00 2022-10-24 09:00:00 1666602000
## 65 0.4097819 2022-10-24 15:00:00 2022-10-24 10:00:00 1666605600
## 66 0.4060566 2022-10-24 16:00:00 2022-10-24 11:00:00 1666609200
## 67 0.3986061 2022-10-24 17:00:00 2022-10-24 12:00:00 1666612800
## 68 0.3911555 2022-10-24 18:00:00 2022-10-24 13:00:00 1666616400
## 69 0.3837049 2022-10-24 19:00:00 2022-10-24 14:00:00 1666620000
## 70 0.3762543 2022-10-24 20:00:00 2022-10-24 15:00:00 1666623600
## 71 0.3725290 2022-10-24 21:00:00 2022-10-24 16:00:00 1666627200
## 72 0.3688037 2022-10-24 22:00:00 2022-10-24 17:00:00 1666630800str(data)## 'data.frame': 72 obs. of 12 variables:
## $ V1 : int 0 1 2 3 4 5 6 7 8 9 ...
## $ AQI : num 30.2 28.2 26.6 25 26 26 26 24 23 24 ...
## $ CO : num 199 198 199 202 205 ...
## $ datetime : chr "2022-10-21:18" "2022-10-21:19" "2022-10-21:20" "2022-10-21:21" ...
## $ no2 : num 0.0469 0.0465 0.0469 0.0482 0.0489 ...
## $ o3 : num 55.8 54.9 54.6 55.1 55.8 ...
## $ pm10 : num 10.5 10.7 11.2 11.1 10.4 ...
## $ pm25 : num 5.64 4.62 3.52 2.23 1.98 ...
## $ so2 : num 0.387 0.41 0.402 0.376 0.339 ...
## $ timestamp_local: POSIXct, format: "2022-10-21 23:00:00" "2022-10-22 00:00:00" ...
## $ timestamp_utc : POSIXct, format: "2022-10-21 18:00:00" "2022-10-21 19:00:00" ...
## $ ts : int 1666375200 1666378800 1666382400 1666386000 1666389600 1666393200 1666396800 1666400400 1666404000 1666407600 ...Classic Regression
regresiklasik = lm(AQI~CO+no2+o3+pm10+pm25+so2, data=data)
summary(regresiklasik)##
## Call:
## lm(formula = AQI ~ CO + no2 + o3 + pm10 + pm25 + so2, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.73339 -0.25673 0.03901 0.24063 0.52176
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.46525 3.49796 -0.419 0.67668
## CO 0.01142 0.01674 0.682 0.49769
## no2 16.32158 5.28647 3.087 0.00297 **
## o3 0.45725 0.01196 38.241 < 2e-16 ***
## pm10 -0.09697 0.02868 -3.381 0.00123 **
## pm25 1.36334 0.08151 16.726 < 2e-16 ***
## so2 -9.45886 1.68392 -5.617 4.36e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3027 on 65 degrees of freedom
## Multiple R-squared: 0.991, Adjusted R-squared: 0.9902
## F-statistic: 1194 on 6 and 65 DF, p-value: < 2.2e-16Dapat terlihat bahwa peubah penjelas yang signifikan terhadap peubah respon adalah x2, x3, x4, x5, dan x6.
rsqklasik = summary(regresiklasik)$r.squared;rsqklasik## [1] 0.9910094rseklasik = sigma(regresiklasik);rseklasik## [1] 0.3027118R-Square yang diperoleh dari data tersebut adalah 0.9910094, artinya peubah-peubah penjelas yang digunakan mampu menjelaskan AQI sebagai peubah respon sebesar 99.10094%. Selain itu, diperoleh pula RSE dari regresi klasik ini sebesar 0.3027118.
Ridge Regression
Peubah yang digunakan
y <- data$AQI
x <- data.matrix(data[,c('CO', 'no2', 'o3', 'pm10', 'pm25', 'so2')])CV
cvr <- cv.glmnet(x,y, alpha=0)
plot(cvr)
Best Model
bestr <- cvr$lambda.min
bestridge <- glmnet(x,y, alpha=0, lambda=bestr)
coef(bestridge)## 7 x 1 sparse Matrix of class "dgCMatrix"
## s0
## (Intercept) -9.25162966
## CO 0.08367469
## no2 34.55907448
## o3 0.33780048
## pm10 -0.05536372
## pm25 1.07011621
## so2 -10.52153546Fungsi R-Square
rsquare <- function(bestmodel,bestlambda,x,y){
yduga <- predict(bestmodel, s=bestlambda, newx=x)
jkt <- sum((y-mean(y))^2)
jkg <- sum((yduga-y)^2)
rsquare <- 1-jkg/jkt
return(rsquare)
}R-Square - Ridge Regression
rsquare(bestridge,bestr,x,y)## [1] 0.9744175R-Square yang diperoleh dari regresi Ridge adalah sebesar 0.9744175. Hasil ini lebih kecil dari R-Square di regresi klasik. Namun, akan dibandingkan pula dengan R-Square regresi Lasso.
RSE - Ridge Regression
# Prediksi model ridge pada data pelatihan
train_predictionsridge <- predict(bestridge,newx = x)
# Hitung residu (selisih antara prediksi dan nilai sebenarnya)
residualsridge <- y - train_predictionsridge
# Hitung varian residu
dfridge <- length(y) - length(bestridge$beta)
residual_varianceridge <- sum(residualsridge^2) / dfridge
# Hitung RSE
rseridge <- sqrt(residual_varianceridge)
# Tampilkan hasil RSE
print(paste("Residual Standard Error (RSE):",rseridge))## [1] "Residual Standard Error (RSE): 0.506746053702993"Diperoleh pula nilai RSE pada regresi Ridge sebesar 0.506746053702993
Lasso Regression
Data Peubah
y <- data$AQI
x <- data.matrix(data[,c('CO', 'no2', 'o3', 'pm10', 'pm25', 'so2')])CV
cvl <- cv.glmnet(x,y,alpha=1)
plot(cvl)
Best Model
bestl <-cvl$lambda.min
bestlasso <- glmnet(x,y,alpha=1,lambda=bestl)
coef(bestlasso)## 7 x 1 sparse Matrix of class "dgCMatrix"
## s0
## (Intercept) -1.79850521
## CO 0.01216225
## no2 15.29552875
## o3 0.45776935
## pm10 -0.08576635
## pm25 1.32314730
## so2 -8.94032751Fungsi R-Square
rsquare <- function(bestmodel,bestlambda,x,y){
yduga <- predict(bestmodel, s=bestlambda, newx=x)
jkt <- sum((y-mean(y))^2)
jkg <- sum((yduga-y)^2)
rsquare <- 1-jkg/jkt
return(rsquare)
}R-Square - Lasso Regression
rsquare(bestlasso,bestl,x,y)## [1] 0.9909643R-Square yang diperoleh dari regresi Lasso adalah 0.9909643. Nilai ini hampir sama dengan R-Square regresi klasik dan lebih besar dari nilai R-Square regresi Ridge.
RSE - Lasso Regression
# Prediksi model Lasso pada data pelatihan
train_predictionslasso <- predict(bestlasso,newx = x)
# Hitung residu (selisih antara prediksi dan nilai sebenarnya)
residualslasso <- y - train_predictionslasso
# Hitung varian residu
dflasso <- length(y) - length(bestlasso$beta)
residual_variancelasso <- sum(residualslasso^2) / dflasso
# Hitung RSE
rselasso <- sqrt(residual_variancelasso)
# Tampilkan hasil RSE
print(paste("Residual Standard Error (RSE):",rselasso))## [1] "Residual Standard Error (RSE): 0.301162226604557"RSE yang diperoleh dari regresi Lasso adalah 0.301162226604557, nilai ini lebih kecil dari RSE regresi Ridge sebelumnya.
Keputusan
Untuk data AQI yang digunakan di sini, regresi Lasso dapat dikatakan yang terbaik, karena walaupun regresi klasik memiliki nilai R-Square yang lebih tinggi (0.991) daripada regresi Lasso (0.9909643), nilai RSE regresi klasik (0.3027) sedikit lebih tinggi dari RSE regresi Lasso (0.301162226604557).
Model Terbaik
bestl <-cvl$lambda.min
bestlasso <- glmnet(x,y,alpha=1,lambda=bestl)
coef(bestlasso)## 7 x 1 sparse Matrix of class "dgCMatrix"
## s0
## (Intercept) -1.79850521
## CO 0.01216225
## no2 15.29552875
## o3 0.45776935
## pm10 -0.08576635
## pm25 1.32314730
## so2 -8.94032751Kesimpulan
Semua peubah (x1, x2, x3, x4, x5, dan x6) yang digunakan bersignifikan terhadap peubah respon (AQI). x1 (CO), x2 (no2), x3 (o3), dan x5 (pm25) berpengaruh positif terhadap peubah respon, artinya semakin tinggi nilai peubah-peubah penjelas tersebut, maka semakin tinggi pula nilai y. Sementara itu, x4 (pm10) dan x6 (so2) berpengaruh negatif terhadap nilai y, artinya semakin tinggi nilai kedua peubah tersebut, maka akan semakin rendah nilai peubah y.