Model Building

Dari bagian sebelumnya, mendapat ilham tentang model dan model building. Prinsip dasar dari membangun sebuah model melibat serangkaian proses tidak berurutan dan kontinuu terkait : * Eksplorasi — > Menemukan pattern * Pemilihan family model (berdasarkan kebutuhan dan kondisi riil data) * Model fitting * Repeat Satu hal yang perlu kita pahami adalah, carilah model paling optimal; yaitu model yang menghasilkan residual paling minimal. Membangun model hampir mirip dengan melukis. Kita harus tahu kapan untuk berhenti mencari kesempurnaan dari karya kita, dan berani menguji cobanya ke lapangan.

Untuk pembahasan ini, kita akan mulai dengan load seluruh data dan package yang kita butuhkan.

Diamonds

Jika kita ingat kembali pada dataset ‘diamonds’, ada satu fakta menarik di mana berlian berkualitas rendah (poor cuts, bad color, dan inferior clarity) memiliki harga lebih mahal daripada berlian berkualitas tinggi.

diamonds %>% ggplot(aes(cut, price)) +geom_boxplot()

diamonds %>% ggplot(aes(color,price)) + geom_boxplot()

diamonds %>% ggplot(aes(clarity, price)) + geom_boxplot()

Kenapa bisa seperti itu?

Jawabannya ada pada eksplorasi data carat vs price.

ggplot(diamonds, aes(carat, price)) + geom_hex(bins = 50)

# Jika bingung, silahkan lihat pada grafik garis berikut
ggplot(diamonds, aes(carat, price)) + geom_line()

Yup, grafik di atas menunjukkan hubungan linier positif antara carat dengan price diamond; dan hubungan itu sangat kuat. Kita mungkin bisa membuat hipotesis bahwa karat merupakan faktor siginigikan dalam penentuan harga diamond.

Namun untuk memperkuat analisa. Mari kita analisa bagaimana kondisi karat dalam berbagai kategori kualitas diamond

ggplot(diamonds, aes(cut, carat)) + geom_boxplot()

ggplot(diamonds, aes(clarity, carat)) + geom_boxplot()

ggplot(diamonds, aes(color, carat)) + geom_boxplot()

Pada boxplot ditunjukkan bahwa diamond dengan carat tinggi terdistribusi lebih banyak pada diamond dengan kualitas renda. Ini memperkuat hipotesis kita bahwa carat merupakan penentu utama dari harga diamond.

Setelah memiliki hipotesis secara visual, mari kita coba buktikan dugaan kita melalui modelling. Kita tidak bisa terlalu cepat mengambil kesimpulan carat sebagai satu-satunya variabel yang berpengaruh terhadap harga bukan?

Nah, sebelum lebih jauh, satu hal yg perlu kita pahami saat berhadapan dengan big data seperti diamonds, ada baiknya kita memperhatikan distribusi masing-masing data. (sebenarnya, lebih sempurna jika kita melakukan ini pada setiap proses analisa data)

Untuk melakukannya, kita bisa menggunakan fungsi summary().

summary(diamonds)
##      carat               cut        color        clarity          depth      
##  Min.   :0.2000   Fair     : 1610   D: 6775   SI1    :13065   Min.   :43.00  
##  1st Qu.:0.4000   Good     : 4906   E: 9797   VS2    :12258   1st Qu.:61.00  
##  Median :0.7000   Very Good:12082   F: 9542   SI2    : 9194   Median :61.80  
##  Mean   :0.7979   Premium  :13791   G:11292   VS1    : 8171   Mean   :61.75  
##  3rd Qu.:1.0400   Ideal    :21551   H: 8304   VVS2   : 5066   3rd Qu.:62.50  
##  Max.   :5.0100                     I: 5422   VVS1   : 3655   Max.   :79.00  
##                                     J: 2808   (Other): 2531                  
##      table           price             x                y         
##  Min.   :43.00   Min.   :  326   Min.   : 0.000   Min.   : 0.000  
##  1st Qu.:56.00   1st Qu.:  950   1st Qu.: 4.710   1st Qu.: 4.720  
##  Median :57.00   Median : 2401   Median : 5.700   Median : 5.710  
##  Mean   :57.46   Mean   : 3933   Mean   : 5.731   Mean   : 5.735  
##  3rd Qu.:59.00   3rd Qu.: 5324   3rd Qu.: 6.540   3rd Qu.: 6.540  
##  Max.   :95.00   Max.   :18823   Max.   :10.740   Max.   :58.900  
##                                                                   
##        z         
##  Min.   : 0.000  
##  1st Qu.: 2.910  
##  Median : 3.530  
##  Mean   : 3.539  
##  3rd Qu.: 4.040  
##  Max.   :31.800  
## 

Dari summary ini kita bisa mengetahu gambaran umum tentang sebaran data. Namun ada satu parameter yang kurang di sini yaitu modus (mode). Untuk itu, mari kita buat terlebih dahulu fungsi modus kita yang kita beri nama modes().

# mode fucntion
modes <- function(x){
  ux <- unique(x)
  tab <- tabulate(match(x,ux))
  ux[tab == max(tab)]
}

Setelah itu, kita mulai lakukan analisa modus.

modes(diamonds$carat)
## [1] 0.3
modes(diamonds$depth)
## [1] 62
modes(diamonds$table)
## [1] 56
modes(diamonds$price)
## [1] 605
modes(diamonds$x)
## [1] 4.37
modes(diamonds$y)
## [1] 4.34
modes(diamonds$z)
## [1] 2.7

Di sini kita bisa mengetahui pattern terkuat berada dari variabel ‘carat’ terhadap ‘price’.

ggplot(diamonds,aes(carat, price)) +
  geom_hex(bins = 50)

Jika kita ingin mengevaluasi dampak attribut diamond lain terhadap harga, kita bisa membangun model untuk memisahkan efek dari ‘carat’ Caranya : * Pertama, kita transformasikan carat dan price menjadi log.

diamonds2 <- diamonds %>% 
  mutate(lprice = log2(price), lcarat = log2(carat))

diamonds2
## # A tibble: 53,940 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.29 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,930 more rows

Apa gunanya transformasi log? Hal ini untuk mempermudah kita melihat hubungan antara carat dan price. Transformasi log akan membuat pattern linier dan memudahkan kita untuk bekerja.

ggplot(diamonds2, aes(lcarat, lprice)) + geom_hex(bins = 50)

Lihat? hasil transformasi log membuat hubungan linier kedua variabel di atas semakin jelas. Langkah selanjutnya dalam mengeluarkan pola linier dari data adalah dengan membuat hubungan carat vs price eksplisit melalui model fitting.

mod_diamonds <- lm(lprice~lcarat, data = diamonds2)
mod_diamonds
## 
## Call:
## lm(formula = lprice ~ lcarat, data = diamonds2)
## 
## Coefficients:
## (Intercept)       lcarat  
##      12.189        1.676
summary(mod_diamonds)
## 
## Call:
## lm(formula = lprice ~ lcarat, data = diamonds2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.17606 -0.24456 -0.00853  0.24002  1.93023 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 12.188841   0.001969  6190.9   <2e-16 ***
## lcarat       1.675817   0.001934   866.6   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3789 on 53938 degrees of freedom
## Multiple R-squared:  0.933,  Adjusted R-squared:  0.933 
## F-statistic: 7.51e+05 on 1 and 53938 DF,  p-value: < 2.2e-16

Langkah selanjutnya adalah dengan mencoba memvisualisasikan model dengan sedikit perubahan.

# Membuat objek baru untuk visualisasi model menggunakan data_grid()
grid <- diamonds2 %>% 
  data_grid(carat) %>% 
  mutate(lcarat= log2(carat)) %>%
  add_predictions(mod_diamonds, "lprice") %>%
  mutate(price = 2 ^ lprice)

# lihat hasilnya
grid
## # A tibble: 273 x 4
##    carat lcarat lprice price
##    <dbl>  <dbl>  <dbl> <dbl>
##  1  0.2   -2.32   8.30  315.
##  2  0.21  -2.25   8.42  341.
##  3  0.22  -2.18   8.53  369.
##  4  0.23  -2.12   8.64  398.
##  5  0.24  -2.06   8.74  427.
##  6  0.25  -2      8.84  457.
##  7  0.26  -1.94   8.93  488.
##  8  0.27  -1.89   9.02  520.
##  9  0.28  -1.84   9.11  553.
## 10  0.29  -1.79   9.20  587.
## # ... with 263 more rows

Sekarang, kita visualisasikan model di atas

ggplot(diamonds2, aes(carat, price)) + 
  geom_hex(bins = 50) +
  geom_line(data = grid, color ="red", size = 1)

Perhatikan bahwa dalam model di atas, kita tidak lagi menggunakan prediksi linier pada variabel ‘price’ (karena setelah mendapat nilai prediksi lprice, kita mengubahnya lagi menjadi ‘price’)

Apakah transformasi ini berdampak ? Kita bisa cek pada residualnya.

# Cek residual
diamonds2 <- diamonds2 %>% add_residuals(mod_diamonds, "lresid")
diamonds2 <- diamonds2 %>% mutate(resid = 2^lresid)
diamonds2
## # A tibble: 53,940 x 14
##    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.29 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,930 more rows, and 2 more variables: lresid <dbl>, resid <dbl>

Untuk lebih mudah, kita lihat visualnya

ggplot(diamonds2, aes(lcarat, resid)) + geom_hex(bins = 50)

ggplot(diamonds2, aes(lcarat, lresid)) + geom_hex(bins = 50)

Kita bisa menggunakan ‘lresid’ pada attribute diamonds lain

ggplot(diamonds2, aes(cut, lresid)) + geom_boxplot()