Rmodel
Arwan Zhagi
2/9/2021
Model Basics with rmodel
Memuat Package yang digunakan
library(tidyverse)
library(modelr)
options(na.action = na.warn)
library(nycflights13)
library(lubridate)Model Building
Dataset yang digunakan
diamonds## # A tibble: 53,940 x 10
## carat cut color clarity depth table price x y z
## <dbl> <ord> <ord> <ord> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 0.23 Ideal E SI2 61.5 55 326 3.95 3.98 2.43
## 2 0.21 Premium E SI1 59.8 61 326 3.89 3.84 2.31
## 3 0.23 Good E VS1 56.9 65 327 4.05 4.07 2.31
## 4 0.290 Premium I VS2 62.4 58 334 4.2 4.23 2.63
## 5 0.31 Good J SI2 63.3 58 335 4.34 4.35 2.75
## 6 0.24 Very Good J VVS2 62.8 57 336 3.94 3.96 2.48
## 7 0.24 Very Good I VVS1 62.3 57 336 3.95 3.98 2.47
## 8 0.26 Very Good H SI1 61.9 55 337 4.07 4.11 2.53
## 9 0.22 Fair E VS2 65.1 61 337 3.87 3.78 2.49
## 10 0.23 Very Good H VS1 59.4 61 338 4 4.05 2.39
## # ... with 53,930 more rowsHubungan antara kualitas dan harga diamond pada dataset diamond.
ggplot(diamonds, aes(cut, price)) +
geom_boxplot()
ggplot(diamonds, aes(color, price)) +
geom_boxplot()
ggplot(diamonds, aes(clarity, price)) +
geom_boxplot()
Diketahui, potongan paling jelek adalah fair, warna paling jelek adalah J, dan clarity paling jelel adalah I1.
Pada boxplot diatas dapat terlihat bahwa seakan-akan harga dari diamonds dengan potongan fair, warna J, dan clarity I1 cukup tinggi. Hal ini terjadi karena tidak dilibatkannya variabel yang lebih superior dalam mempengaruhi nilai harga, yaitu variabel carat/weight.
ggplot(diamonds, aes(carat, price)) +
geom_hex(bins = 50)
Mari lihat bagaimana variabel lain dari berlian dapat memengruhi harga dengan menyesuaikan model untuk mengurangi efek karat.
Dataset diambil atau difilter dengan carat dibawah 2.5, karena sebagian besar data berada pada karat dibawah 2.5. Hal ini dapat memudahkan.
Selain itu pada visualisasi dataset terlihat bahwa pertambahan carat mengakibatkan harga naik pesat (seperti grafik eksponen). Supaya grafik menjadi linear lakukan transformasi log.
diamonds2 <- diamonds %>%
filter(carat <= 2.5) %>%
mutate(lprice = log2(price),
lcarat = log2(carat))
diamonds2## # A tibble: 53,814 x 12
## carat cut color clarity depth table price x y z lprice lcarat
## <dbl> <ord> <ord> <ord> <dbl> <dbl> <int> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0.23 Ideal E SI2 61.5 55 326 3.95 3.98 2.43 8.35 -2.12
## 2 0.21 Premium E SI1 59.8 61 326 3.89 3.84 2.31 8.35 -2.25
## 3 0.23 Good E VS1 56.9 65 327 4.05 4.07 2.31 8.35 -2.12
## 4 0.290 Premium I VS2 62.4 58 334 4.2 4.23 2.63 8.38 -1.79
## 5 0.31 Good J SI2 63.3 58 335 4.34 4.35 2.75 8.39 -1.69
## 6 0.24 Very G~ J VVS2 62.8 57 336 3.94 3.96 2.48 8.39 -2.06
## 7 0.24 Very G~ I VVS1 62.3 57 336 3.95 3.98 2.47 8.39 -2.06
## 8 0.26 Very G~ H SI1 61.9 55 337 4.07 4.11 2.53 8.40 -1.94
## 9 0.22 Fair E VS2 65.1 61 337 3.87 3.78 2.49 8.40 -2.18
## 10 0.23 Very G~ H VS1 59.4 61 338 4 4.05 2.39 8.40 -2.12
## # ... with 53,804 more rowsVisualisasi dari harga dan carat setelah ditasformasi.
ggplot(diamonds2, aes(lcarat, lprice)) +
geom_hex(bins = 50)
Model yang terbentuk dari transformasi harga dan carat.
mod_diamond <- lm(lprice ~ lcarat, data = diamonds2)
mod_diamond##
## Call:
## lm(formula = lprice ~ lcarat, data = diamonds2)
##
## Coefficients:
## (Intercept) lcarat
## 12.194 1.681Hasil prediksi harga dari 50 data berdasarkan model yang dibentuk.
grid <- diamonds2 %>%
data_grid(carat = seq_range(carat, 50)) %>%
mutate(lcarat = log2(carat)) %>%
add_predictions(mod_diamond, "lprice") %>%
mutate(price = 2 ^ lprice)
grid## # A tibble: 50 x 4
## carat lcarat lprice price
## <dbl> <dbl> <dbl> <dbl>
## 1 0.2 -2.32 8.29 313.
## 2 0.247 -2.02 8.80 446.
## 3 0.294 -1.77 9.22 598.
## 4 0.341 -1.55 9.58 767.
## 5 0.388 -1.37 9.90 953.
## 6 0.435 -1.20 10.2 1154.
## 7 0.482 -1.05 10.4 1372.
## 8 0.529 -0.920 10.6 1604.
## 9 0.576 -0.797 10.9 1850.
## 10 0.622 -0.684 11.0 2111.
## # ... with 40 more rowsVisualisasi model dengan dataset.
ggplot(diamonds2, aes(carat, price)) +
geom_hex(bins = 50) +
geom_line(data = grid, color = "red", size = 2)
Visualisasi dari residual. Terlihat bahwa kecenderungan linearitas dari residual sudah dihilangkan. Pada visualisasi residual ini juga didapat bahwa kebanyakan diamond lebih murah dari hasil dari model (expected).
diamonds2 <- diamonds2 %>%
add_residuals(mod_diamond, "lresid")
ggplot(diamonds2, aes(lcarat, lresid)) +
geom_hex(bins = 50)
sebagai perbandingan berikut visualisasi dari residual tanpa di transformasi.
mod_diamond1 <- lm(price ~ carat, data = diamonds2)
grid1 <- diamonds2 %>%
data_grid(carat = seq_range(carat, 50)) %>%
add_predictions(mod_diamond1, "price")
diamonds1 <- diamonds2 %>%
add_residuals(mod_diamond1, "resid")
ggplot(diamonds1, aes(carat, resid)) +
geom_hex(bins = 50)
Berdasarkan residual tersebut dapat dilihat perbandingan harga diamond dari tiap masing-masing cut, color dan clarity.
ggplot(diamonds2, aes(cut, lresid)) +
geom_boxplot()
ggplot(diamonds2, aes(color, lresid)) +
geom_boxplot() 
ggplot(diamonds2, aes(clarity, lresid)) +
geom_boxplot()
A More Complicated Model
memasukkan variabel Cut, Color, dan Clarity pada model. Model yang diguakan model yang eksplisit.
** Masih belum ngerti sih kenapa langsung menggunakan plus, terus si crossnya gak diperiksa.
Model.
mod_diamond2 <- lm( lprice ~ lcarat + color + cut + clarity, data = diamonds2 )
mod_diamond2##
## Call:
## lm(formula = lprice ~ lcarat + color + cut + clarity, data = diamonds2)
##
## Coefficients:
## (Intercept) lcarat color.L color.Q color.C color^4
## 12.206978 1.886239 -0.633998 -0.137580 -0.022072 0.016570
## color^5 color^6 cut.L cut.Q cut.C cut^4
## -0.002828 0.003533 0.173866 -0.050346 0.019129 -0.002410
## clarity.L clarity.Q clarity.C clarity^4 clarity^5 clarity^6
## 1.308155 -0.334090 0.178423 -0.088059 0.035885 -0.001371
## clarity^7
## 0.048221Karena predictornya banyak, jadi susah kalau langsung divisualisasikan. Jadi dipisah per variabel bebas.
Prediktor Cut.
Predictor Cut diambil semua valuenya, sedangkan preictor lain apabila prediktor tersebut kontinu diambil mediannya, apabila prediktor tersebut kategorik diambil modusnya.
grid <- diamonds2 %>%
data_grid(cut, .model = mod_diamond2) %>%
add_predictions(mod_diamond2)
grid ## # A tibble: 5 x 5
## cut lcarat color clarity pred
## <ord> <dbl> <chr> <chr> <dbl>
## 1 Fair -0.515 G VS2 11.2
## 2 Good -0.515 G VS2 11.3
## 3 Very Good -0.515 G VS2 11.4
## 4 Premium -0.515 G VS2 11.4
## 5 Ideal -0.515 G VS2 11.4Visualisasi dari prediktor Cut menggunakan model dengan kemepat prediktor.
ggplot(grid, aes(cut, pred)) +
geom_point(color="red", size=4)
Predictor Color
grid <- diamonds2 %>%
data_grid(color, .model = mod_diamond2) %>%
add_predictions(mod_diamond2)
grid ## # A tibble: 7 x 5
## color lcarat cut clarity pred
## <ord> <dbl> <chr> <chr> <dbl>
## 1 D -0.515 Premium VS2 11.6
## 2 E -0.515 Premium VS2 11.6
## 3 F -0.515 Premium VS2 11.5
## 4 G -0.515 Premium VS2 11.4
## 5 H -0.515 Premium VS2 11.3
## 6 I -0.515 Premium VS2 11.1
## 7 J -0.515 Premium VS2 10.9Visualisasi dari prediktor Color menggunakan model dengan kemepat prediktor.
ggplot(grid, aes(color, pred)) +
geom_point(color="red", size=4)
Predictor Clarity
grid <- diamonds2 %>%
data_grid(clarity, .model = mod_diamond2) %>%
add_predictions(mod_diamond2)
grid ## # A tibble: 8 x 5
## clarity lcarat color cut pred
## <ord> <dbl> <chr> <chr> <dbl>
## 1 I1 -0.515 G Premium 10.4
## 2 SI2 -0.515 G Premium 11.0
## 3 SI1 -0.515 G Premium 11.2
## 4 VS2 -0.515 G Premium 11.4
## 5 VS1 -0.515 G Premium 11.5
## 6 VVS2 -0.515 G Premium 11.7
## 7 VVS1 -0.515 G Premium 11.8
## 8 IF -0.515 G Premium 11.9Visualisasi dari prediktor Clarity menggunakan model dengan kemepat prediktor.
ggplot(grid, aes(clarity, pred)) +
geom_point(color="red", size=4)
RESIDUAL
diamonds2 <- diamonds2 %>%
add_residuals(mod_diamond2, "lresid2")
ggplot(diamonds2, aes(lcarat, lresid2)) +
geom_hex(bins = 50)
Terlihat bahwa terdapat data dengan residual lebih dari 1. misal, pada visualisasi tersebut terdapat data dengan residual 2. Hal ini berarti bahwa harga berlian 4 kali lipat dari harga prediksi (expected)
diamonds2 %>%
filter(abs(lresid2) > 1) %>%
add_predictions(mod_diamond2) %>%
mutate(pred = round(2 ^ pred)) %>%
select(price, pred, carat:table, x:z) %>%
arrange(price) ## # A tibble: 16 x 11
## price pred carat cut color clarity depth table x y z
## <int> <dbl> <dbl> <ord> <ord> <ord> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1013 264 0.25 Fair F SI2 54.4 64 4.3 4.23 2.32
## 2 1186 284 0.25 Premium G SI2 59 60 5.33 5.28 3.12
## 3 1186 284 0.25 Premium G SI2 58.8 60 5.33 5.28 3.12
## 4 1262 2644 1.03 Fair E I1 78.2 54 5.72 5.59 4.42
## 5 1415 639 0.35 Fair G VS2 65.9 54 5.57 5.53 3.66
## 6 1415 639 0.35 Fair G VS2 65.9 54 5.57 5.53 3.66
## 7 1715 576 0.32 Fair F VS2 59.6 60 4.42 4.34 2.61
## 8 1776 412 0.290 Fair F SI1 55.8 60 4.48 4.41 2.48
## 9 2160 314 0.34 Fair F I1 55.8 62 4.72 4.6 2.6
## 10 2366 774 0.3 Very Good D VVS2 60.6 58 4.33 4.35 2.63
## 11 3360 1373 0.51 Premium F SI1 62.7 62 5.09 4.96 3.15
## 12 3807 1540 0.61 Good F SI2 62.5 65 5.36 5.29 3.33
## 13 3920 1705 0.51 Fair F VVS2 65.4 60 4.98 4.9 3.23
## 14 4368 1705 0.51 Fair F VVS2 60.7 66 5.21 5.11 3.13
## 15 10011 4048 1.01 Fair D SI2 64.6 58 6.25 6.2 4.02
## 16 10470 23622 2.46 Premium E SI2 59.7 59 8.82 8.76 5.25MENCOBA MENGGUNKAN DATA DENGAN MEMBUANG DATA YANG DIATAS
diamonds2<-diamonds2 %>%
filter(abs(lresid2) <= 1) Model.
mod_diamond2 <- lm( lprice ~ lcarat + color + cut + clarity, data = diamonds2 )
mod_diamond2##
## Call:
## lm(formula = lprice ~ lcarat + color + cut + clarity, data = diamonds2)
##
## Coefficients:
## (Intercept) lcarat color.L color.Q color.C color^4
## 12.205752 1.887263 -0.634200 -0.137159 -0.022361 0.015946
## color^5 color^6 cut.L cut.Q cut.C cut^4
## -0.001845 0.003202 0.179443 -0.054694 0.021984 -0.003214
## clarity.L clarity.Q clarity.C clarity^4 clarity^5 clarity^6
## 1.309477 -0.334049 0.178549 -0.087779 0.035417 -0.001623
## clarity^7
## 0.047751Prediktor Cut.
grid <- diamonds2 %>%
data_grid(cut, .model = mod_diamond2) %>%
add_predictions(mod_diamond2)
grid ## # A tibble: 5 x 5
## cut lcarat color clarity pred
## <ord> <dbl> <chr> <chr> <dbl>
## 1 Fair -0.515 G VS2 11.2
## 2 Good -0.515 G VS2 11.3
## 3 Very Good -0.515 G VS2 11.4
## 4 Premium -0.515 G VS2 11.4
## 5 Ideal -0.515 G VS2 11.4Visualisasi dari prediktor Cut menggunakan model dengan kemepat prediktor.
ggplot(grid, aes(cut, pred)) +
geom_point(color="red", size=4)
Predictor Color
grid <- diamonds2 %>%
data_grid(color, .model = mod_diamond2) %>%
add_predictions(mod_diamond2)
grid ## # A tibble: 7 x 5
## color lcarat cut clarity pred
## <ord> <dbl> <chr> <chr> <dbl>
## 1 D -0.515 Premium VS2 11.6
## 2 E -0.515 Premium VS2 11.6
## 3 F -0.515 Premium VS2 11.5
## 4 G -0.515 Premium VS2 11.4
## 5 H -0.515 Premium VS2 11.3
## 6 I -0.515 Premium VS2 11.1
## 7 J -0.515 Premium VS2 10.9Visualisasi dari prediktor Color menggunakan model dengan kemepat prediktor.
ggplot(grid, aes(color, pred)) +
geom_point(color="red", size=4)
Predictor Clarity
grid <- diamonds2 %>%
data_grid(clarity, .model = mod_diamond2) %>%
add_predictions(mod_diamond2)
grid ## # A tibble: 8 x 5
## clarity lcarat color cut pred
## <ord> <dbl> <chr> <chr> <dbl>
## 1 I1 -0.515 G Premium 10.4
## 2 SI2 -0.515 G Premium 11.0
## 3 SI1 -0.515 G Premium 11.2
## 4 VS2 -0.515 G Premium 11.4
## 5 VS1 -0.515 G Premium 11.5
## 6 VVS2 -0.515 G Premium 11.7
## 7 VVS1 -0.515 G Premium 11.8
## 8 IF -0.515 G Premium 11.9Visualisasi dari prediktor Clarity menggunakan model dengan kemepat prediktor.
ggplot(grid, aes(clarity, pred)) +
geom_point(color="red", size=4)
RESIDUAL
diamonds2 <- diamonds2 %>%
add_residuals(mod_diamond2, "lresid2")
ggplot(diamonds2, aes(lcarat, lresid2)) +
geom_hex(bins = 50)
##SAMA AJAH TIDAK TERLALU BERARTIā¦
##BERSAMBUNG