Modelling dalam statistik dan matematika adalah proses belajar yang panjang dan sangat mungkin membingungkan. Selama ini, kita mengenal modelling sebagai upaya sistematis untuk menjabarkan hubungan matematis-statistik antar satu dengan beberapa variabel, dan menghasilkan suatu output tertentu (biasanya berupa prediksi atau sinyal tertentu).

Banyak terminologi yang bertebaran dari berbagai sumber. Namun secara umum, terdapat dua jenis permodelan : * Data-finding model (unsupervised model) yang berguna untuk mencari pola dan hubungan dalam model. * Predictive model (supervised model) yang berguna untuk melakukan prediksi atas serangkaian variabel yang berhubungan.

Secara general, model bertujuan untuk melakukan inferensi, konfirmasi hipotesis, atau eksplorasi.

Untuk melakukan analisa konfirmasi hipotesis, salah satu pendekatan yang bisa dipakai adalah dengan membagi dataset kita menjadi 3 bagian : * 60% untuk training data (explorasi) di mana kita bisa bebas menggunakan bagian data ini untuk mengeksplorasi data * 20% untuk query di mana kita menggunakan untuk mengomparasikan model-model. * 20% untuk test; di mana bagian ini hanya bisa digunakan sekali.

Dalam pembentukan sebuah model, kita akan selalu mengenal yang namanya Family model dan fitted model. Family model adalah struktur persis tapi generik yang mampu memberikan gambaran pola yang di cari. Contoh : Pola garis lurus dari data bisa jadi menunjukkan bahwa data yang ada model family linier. Sedangkan fitted model adalah model dari family model yang paling akurat untuk data yang kita miliki.

Seperti yang dijelaskan di atas, dalam modelling, hal pertama yang perlu kita analisa adalah pada bagian model family apa sekiranya data kita termasuk.

Contoh, menggunakan dataset ‘sim1’.

sim1
## # A tibble: 30 x 2
##        x     y
##    <int> <dbl>
##  1     1  4.20
##  2     1  7.51
##  3     1  2.13
##  4     2  8.99
##  5     2 10.2 
##  6     2 11.3 
##  7     3  7.36
##  8     3 10.5 
##  9     3 10.5 
## 10     4 12.4 
## # ... with 20 more rows

Pertama, kita tahu bahwa dataset sim1 terdiri dari 2 variabel x dan y yang sama-sama memiliki data numeric. Sekarang, cara termudah untuk melihat pola antara x dan y adalah dengan memvisualisasikannya.

sim1 %>% ggplot(aes(x,y)) + geom_point()

Pada data di atas, kita sudah bisa melihat atau menerka secara visual termasuk family model apakah variabel x~y. Yep, data di atas menunjukkan adanya relasi linear antara x dan y yang bisa digambarkan dengan family model : y = a_0 + a_1*x.

Tapi a_0 dan a_1 mana yang paling ‘tepat’ untuk menjadi parameter model? Untuk memberikan logika yang pas, kita akan lakukan beberapa simulasi.

models <- tibble(
  a1 = runif(250, -20, 40),
  a2 = runif(250, -5, 5)
)

ggplot(sim1, aes(x,y)) + geom_abline(aes(intercept = a1, slope = a2), data = models, alpha = 1/4) + geom_point()

Graph di atas adalah simulasi intercept (a_0) dan slop (a_1) menggunakan 250 titip intersep dan slop dari objek models. Berbekal siimulasi di atas, kita bisa mencari mana sekiranya garis yang paling mendekati sebaran data.

Logikanya, model fit terbaik dimiliki oleh model yang memiliki varians terkecil. Bahasa mudahnya, mana model yang mampu menghasilkan selisih antara y value dari model (prediction) dengan y value dari data (response) terkecil, itulah model terbaik untuk data tersebut.

Sekedar penajaman logika, kita coba menghitung seluruh nilai beda antara y_model - y_data atau yang kita sebut dengan diff_y. Sebelum itu, kita perlu menghitung seluruh nilai y_model.

# Kita tahu bahwa family modelnya y_model = a_0 + a_1 * x (model linear)
# Kita hitung nilai y_modelnya terlebih dahulu.

model1 <- function(a, data){
  a[1] + a[2]*data$x # Pembuatan fungsi untuk menggambar family model
}

# Kita uji fungsi di atas, dengan menggunakan c(7,1.5)
model1(c(7,1.5),sim1)
##  [1]  8.5  8.5  8.5 10.0 10.0 10.0 11.5 11.5 11.5 13.0 13.0 13.0 14.5 14.5 14.5
## [16] 16.0 16.0 16.0 17.5 17.5 17.5 19.0 19.0 19.0 20.5 20.5 20.5 22.0 22.0 22.0

Kita sudah punya fungsi untuk menghitung nilai y_model. Sekarang, kita perlu membuat fungsi untuk menghitung selisih antara y_model dengan y_data.

# kita bisa menulisnya seperti ini :

measure_dist <- function(mod,data){
diff <- data$y - model1(mod,data)
}

Akan tetapi fungsi di atas akan menghasilkan 30 nilai selisih antara y_model dan y_data. Kita tidak bisa berbuat dengan hanya menggunakan selisih tersebut. Dalam menganalisa mana model terbaik, kita perlu menghitung berapa nilai keseluruhan yang menjadi selisih antara y_model dengan y_data.

Dalam statistik, kita bisa menggunkan konsep root mean squared untuk membuat keseluruhan selisih menjadi satu nilai yang representatif.

Kita ubah fungsi di atas, menjadi seperti ini :

# Fungsi selisih
dist <- function(mod, data){
  diff <- data$y-model1(mod,data)
  diff
}
dist(c(7,1.5),sim1)
##  [1] -4.30008703 -0.98936589 -6.37452722 -1.01114257  0.24310544  1.29682321
##  [7] -4.14363532 -0.99465060 -0.98839921 -0.56541090 -1.10739877  1.25796408
## [13]  4.63004979 -2.76197880  1.52485390 -2.72602298 -0.04402503  0.89479618
## [19]  2.58599269 -0.32814965  2.43630884  2.72590251 -0.60908709  3.47555264
## [25]  6.27700986  2.30510979  0.62830529  2.96809938  1.34642209 -0.02479930
# fungsi standart distance
measure_dist <- function(mod, data){
  diff <- data$y - model1(mod, data)
  sqrt(mean(diff^2))
}

measure_dist(c(7,1.5),sim1)
## [1] 2.665212

Kita sudah punya dua fungsi yang menjadi dasar pembentukan model linear kita. Lalu sekarang bagaimana? Sekarang, kita terapkan seluruh fungsi di atas untuk mencari mana model terbaik yang paling dekat dengan data; dilihat dari selisih terkecilnya.

sim1_dist <- function(a1,a2){
  measure_dist(c(a1,a2), sim1)
}
# Kenapa kita membuat fungsi lagi untuk membungkus fungsi measure_dist? Jawabannya karena fungsi measure_dist hanya mensyaratkan vector dengan panjang 2

#sekarang, kita gunakan fungsi sim1_dist untuk membuat tibble yang berisi kemungkinan intersep dan slop dari a1 dan2 objek 'models'

models <- models %>% 
  mutate(dist = map2_dbl(a1,a2,sim1_dist))
models
## # A tibble: 250 x 3
##       a1     a2  dist
##    <dbl>  <dbl> <dbl>
##  1 34.4   2.53  32.9 
##  2  3.25  2.78   4.23
##  3  6.92 -0.199 11.8 
##  4 -6.24 -2.63  38.7 
##  5 35.3   3.67  40.3 
##  6  4.28  4.11  13.0 
##  7  8.14 -1.59  19.3 
##  8 -5.21 -2.63  37.7 
##  9 28.3   2.32  25.7 
## 10 -6.07  2.50   8.22
## # ... with 240 more rows

Sekarang kita sudah punya 250 calon intersep dan slope untuk model linear kita. Pertanyaannya, bagaimana kita memilih model mana yang terbaik? Jawabannya sederhana, kita tinggal memilih intersep dan slope yang menghasilkan selisih y_model dan y_data paling kecil di antara semuanya. Jika begitu, kita bisa saja langsung memilih model terbaik dengan mengurutkan tibble di atas.

models %>% filter(rank(dist)<=10) %>% arrange(dist)
## # A tibble: 10 x 3
##        a1    a2  dist
##     <dbl> <dbl> <dbl>
##  1  2.34  2.43   2.40
##  2  4.32  1.82   2.50
##  3  6.57  1.56   2.57
##  4  1.65  2.65   2.82
##  5  0.122 2.35   3.38
##  6  1.51  2.92   3.89
##  7 -1.89  2.57   4.17
##  8 11.3   0.780  4.23
##  9  3.25  2.78   4.23
## 10 10.1   0.739  4.53

Pencarian kita bisa berakhir di sini. Tapi merupakan sebuah praktik yang baik (sekaligus membantu menyakinkan klien) jika kita bisa memvisualisasikan 10 model terbaik ini untuk menunjukkan bahwa mereka benar-benar dapat merepresentasikan data yang kita miliki.

Visualisasi 10 model ‘terbaik’ sim1 :

ggplot(sim1, aes(x,y)) +
  geom_point(size = 2, color = "grey30") +
  geom_abline(aes(intercept = a1, slope = a2, color = -dist),
              data = filter(models, rank(dist) <= 10))

Cara lain untuk visualisasi calon-calon model ini adalah membuat visualisasi dari rangkaian intersep dan slope.

ggplot(models, aes(a1,a2)) +
  geom_point(data = filter(models, rank(dist)<=10), size =4, color = "red") +
  geom_point(aes(color = -dist))

Ada cara untuk mengoptimalkan calon-calon model di atas, yaitu dengan menggunakan metode sistematis bernama grid search. Untuk ini, kita perlu membuat calon parameter model baru yang lebih sempit dengan berpatokan pada parameter model yang ditemukan di atas.

# Buat dataframe baru dengan expand.grid()

grid <- expand.grid(
  b1 = seq(-10,15,length = 25),
  b2 = seq(0,4,length = 25)
) %>% mutate(dist = map2_dbl(b1,b2,sim1_dist))
grid %>% arrange(dist)
##              b1        b2      dist
## 1     4.5833333 2.0000000  2.134787
## 2     3.5416667 2.1666667  2.154210
## 3     5.6250000 1.8333333  2.227912
## 4     2.5000000 2.3333333  2.283311
## 5     3.5416667 2.0000000  2.340435
## 6     4.5833333 1.8333333  2.371391
## 7     4.5833333 2.1666667  2.372767
## 8     5.6250000 2.0000000  2.409800
## 9     2.5000000 2.1666667  2.412749
## 10    6.6666667 1.6666667  2.420631
## 11    3.5416667 2.3333333  2.437724
## 12    5.6250000 1.6666667  2.501788
## 13    1.4583333 2.5000000  2.505191
## 14    6.6666667 1.8333333  2.544372
## 15    1.4583333 2.3333333  2.579664
## 16    2.5000000 2.5000000  2.597028
## 17    7.7083333 1.5000000  2.691635
## 18    6.6666667 1.5000000  2.717346
## 19    7.7083333 1.6666667  2.762265
## 20    0.4166667 2.6666667  2.797863
## 21    0.4166667 2.5000000  2.824456
## 22    1.4583333 2.6666667  2.834818
## 23    3.5416667 1.8333333  2.907497
## 24    2.5000000 2.0000000  2.927131
## 25    5.6250000 2.1666667  2.964720
## 26    4.5833333 1.6666667  2.971291
## 27    4.5833333 2.3333333  2.973486
## 28    7.7083333 1.3333333  2.999763
## 29    8.7500000 1.3333333  3.019922
## 30    1.4583333 2.1666667  3.028570
## 31    6.6666667 2.0000000  3.037613
## 32    8.7500000 1.5000000  3.045647
## 33    3.5416667 2.5000000  3.063211
## 34    5.6250000 1.5000000  3.113387
## 35   -0.6250000 2.6666667  3.128902
## 36    0.4166667 2.8333333  3.133275
## 37   -0.6250000 2.8333333  3.141604
## 38    7.7083333 1.8333333  3.186566
## 39    0.4166667 2.3333333  3.204053
## 40    2.5000000 2.6666667  3.227148
## 41    6.6666667 1.3333333  3.323759
## 42    8.7500000 1.1666667  3.332083
## 43    9.7916667 1.3333333  3.378078
## 44    9.7916667 1.1666667  3.388884
## 45    8.7500000 1.6666667  3.401602
## 46   -0.6250000 2.5000000  3.442276
## 47    1.4583333 2.8333333  3.454749
## 48   -0.6250000 3.0000000  3.476812
## 49   -1.6666667 2.8333333  3.477367
## 50   -1.6666667 3.0000000  3.521491
## 51    7.7083333 1.1666667  3.590425
## 52    2.5000000 1.8333333  3.667932
## 53    9.7916667 1.5000000  3.671127
## 54    3.5416667 1.6666667  3.683522
## 55    9.7916667 1.0000000  3.700885
## 56    1.4583333 2.0000000  3.718695
## 57    5.6250000 2.3333333  3.729730
## 58   -1.6666667 2.6666667  3.731240
## 59    0.4166667 3.0000000  3.734391
## 60   10.8333333 1.1666667  3.746524
## 61    6.6666667 2.1666667  3.757558
## 62    4.5833333 1.5000000  3.764642
## 63    4.5833333 2.5000000  3.767241
## 64   10.8333333 1.0000000  3.786652
## 65    0.4166667 2.1666667  3.833178
## 66    7.7083333 2.0000000  3.849310
## 67   -1.6666667 3.1666667  3.853389
## 68   -2.7083333 3.0000000  3.857943
## 69    3.5416667 2.6666667  3.868192
## 70    8.7500000 1.0000000  3.901860
## 71    5.6250000 1.3333333  3.907212
## 72   -2.7083333 3.1666667  3.927048
## 73   10.8333333 1.3333333  3.984098
## 74    8.7500000 1.8333333  4.000589
## 75   -0.6250000 2.3333333  4.005919
## 76    2.5000000 2.8333333  4.027816
## 77   -0.6250000 3.1666667  4.055324
## 78   -2.7083333 2.8333333  4.060126
## 79   10.8333333 0.8333333  4.096329
## 80    6.6666667 1.1666667  4.104834
## 81   11.8750000 1.0000000  4.141383
## 82    9.7916667 1.6666667  4.204976
## 83   11.8750000 0.8333333  4.205057
## 84   -1.6666667 2.5000000  4.229789
## 85    1.4583333 3.0000000  4.239490
## 86    9.7916667 0.8333333  4.248230
## 87   -2.7083333 3.3333333  4.254243
## 88   -3.7500000 3.1666667  4.262035
## 89   11.8750000 1.1666667  4.331106
## 90    7.7083333 1.0000000  4.350013
## 91   -3.7500000 3.3333333  4.351103
## 92   -1.6666667 3.3333333  4.408538
## 93   -3.7500000 3.0000000  4.420032
## 94   10.8333333 1.5000000  4.455168
## 95    0.4166667 3.1666667  4.495868
## 96   -2.7083333 2.6666667  4.497157
## 97   11.8750000 0.6666667  4.511414
## 98    2.5000000 1.6666667  4.525289
## 99    1.4583333 1.8333333  4.541370
## 100  12.9166667 0.8333333  4.555792
## 101   3.5416667 1.5000000  4.563114
## 102   6.6666667 2.3333333  4.598944
## 103   5.6250000 2.5000000  4.601205
## 104   0.4166667 2.0000000  4.610795
## 105  10.8333333 0.6666667  4.621686
## 106   8.7500000 0.8333333  4.635206
## 107  12.9166667 0.6666667  4.638520
## 108   7.7083333 2.1666667  4.649630
## 109   4.5833333 1.3333333  4.653534
## 110   4.5833333 2.6666667  4.656338
## 111  -3.7500000 3.5000000  4.673129
## 112  -4.7916667 3.3333333  4.683562
## 113  12.9166667 1.0000000  4.704626
## 114  -0.6250000 2.1666667  4.731216
## 115  11.8750000 1.3333333  4.743922
## 116   8.7500000 2.0000000  4.751571
## 117   3.5416667 2.8333333  4.762505
## 118  -2.7083333 3.5000000  4.786893
## 119  -4.7916667 3.5000000  4.788744
## 120  -0.6250000 3.3333333  4.789777
## 121   5.6250000 1.1666667  4.793571
## 122  -3.7500000 2.8333333  4.800761
## 123  -4.7916667 3.1666667  4.803990
## 124  -1.6666667 2.3333333  4.898873
## 125   9.7916667 1.8333333  4.901568
## 126   2.5000000 3.0000000  4.916402
## 127  12.9166667 0.5000000  4.941193
## 128   9.7916667 0.6666667  4.953509
## 129   6.6666667 1.0000000  4.979041
## 130  13.9583333 0.6666667  4.984879
## 131  11.8750000 0.5000000  5.016183
## 132  12.9166667 1.1666667  5.064648
## 133  13.9583333 0.5000000  5.083189
## 134  10.8333333 1.6666667  5.095380
## 135  13.9583333 0.8333333  5.098835
## 136  -4.7916667 3.6666667  5.105612
## 137  -2.7083333 2.5000000  5.109119
## 138   1.4583333 3.1666667  5.113721
## 139  -1.6666667 3.5000000  5.114752
## 140  -5.8333333 3.5000000  5.118216
## 141  -4.7916667 3.0000000  5.134178
## 142  -3.7500000 3.6666667  5.184888
## 143   7.7083333 0.8333333  5.205090
## 144  -5.8333333 3.3333333  5.206682
## 145  -5.8333333 3.6666667  5.236566
## 146  10.8333333 0.5000000  5.298957
## 147  11.8750000 1.5000000  5.328227
## 148   0.4166667 3.3333333  5.349661
## 149  -3.7500000 2.6666667  5.356941
## 150  13.9583333 0.3333333  5.382146
## 151  13.9583333 1.0000000  5.411664
## 152  15.0000000 0.5000000  5.425162
## 153  12.9166667 0.3333333  5.427134
## 154   1.4583333 1.6666667  5.436753
## 155   2.5000000 1.5000000  5.444414
## 156  -2.7083333 3.6666667  5.465254
## 157   8.7500000 0.6666667  5.466686
## 158   0.4166667 1.8333333  5.473984
## 159  -5.8333333 3.1666667  5.491979
## 160   3.5416667 1.3333333  5.496779
## 161   6.6666667 2.5000000  5.506379
## 162  15.0000000 0.6666667  5.509293
## 163   7.7083333 2.3333333  5.528094
## 164   5.6250000 2.6666667  5.529031
## 165  15.0000000 0.3333333  5.536365
## 166  -5.8333333 3.8333333  5.548513
## 167  -0.6250000 2.0000000  5.555205
## 168  -6.8750000 3.6666667  5.562923
## 169   4.5833333 1.1666667  5.592594
## 170   8.7500000 2.1666667  5.593659
## 171  12.9166667 1.3333333  5.595239
## 172   4.5833333 2.8333333  5.595511
## 173  -4.7916667 3.8333333  5.598337
## 174  -0.6250000 3.5000000  5.619357
## 175  -6.8750000 3.5000000  5.624086
## 176  -4.7916667 2.8333333  5.637387
## 177  11.8750000 0.3333333  5.666588
## 178  -1.6666667 2.1666667  5.678529
## 179  -6.8750000 3.8333333  5.692167
## 180   9.7916667 2.0000000  5.701561
## 181   3.5416667 3.0000000  5.704286
## 182   5.6250000 1.0000000  5.729679
## 183   9.7916667 0.5000000  5.758987
## 184  15.0000000 0.8333333  5.780236
## 185  15.0000000 0.1666667  5.831740
## 186  -3.7500000 3.8333333  5.836684
## 187  -2.7083333 2.3333333  5.841290
## 188  10.8333333 1.8333333  5.849458
## 189  13.9583333 0.1666667  5.851074
## 190   2.5000000 3.1666667  5.853000
## 191  -6.8750000 3.3333333  5.869708
## 192  13.9583333 1.1666667  5.891772
## 193   6.6666667 0.8333333  5.905160
## 194  -1.6666667 3.6666667  5.918198
## 195  -5.8333333 3.0000000  5.945840
## 196  -6.8750000 4.0000000  5.999525
## 197  -7.9166667 3.8333333  6.015452
## 198  -5.8333333 4.0000000  6.024057
## 199  11.8750000 1.6666667  6.034410
## 200   1.4583333 3.3333333  6.038701
## 201  -3.7500000 2.5000000  6.040300
## 202  12.9166667 0.1666667  6.052362
## 203  -7.9166667 3.6666667  6.053159
## 204  10.8333333 0.3333333  6.077564
## 205   7.7083333 0.6666667  6.115733
## 206  -7.9166667 4.0000000  6.153819
## 207  15.0000000 1.0000000  6.213601
## 208  -4.7916667 4.0000000  6.225300
## 209  -2.7083333 3.8333333  6.241998
## 210  12.9166667 1.5000000  6.253130
## 211   0.4166667 3.5000000  6.258099
## 212  -7.9166667 3.5000000  6.263761
## 213  -4.7916667 2.6666667  6.272110
## 214  -6.8750000 3.1666667  6.278176
## 215  15.0000000 0.0000000  6.285374
## 216   8.7500000 0.5000000  6.357913
## 217   1.4583333 1.5000000  6.374277
## 218   0.4166667 1.6666667  6.388150
## 219   2.5000000 1.3333333  6.398745
## 220  11.8750000 0.1666667  6.418507
## 221  -0.6250000 1.8333333  6.440116
## 222   6.6666667 2.6666667  6.452056
## 223   7.7083333 2.5000000  6.452866
## 224  13.9583333 0.0000000  6.453024
## 225   3.5416667 1.1666667  6.461118
## 226  -8.9583333 4.0000000  6.474165
## 227   5.6250000 2.8333333  6.489081
## 228   8.7500000 2.3333333  6.491495
## 229  -8.9583333 3.8333333  6.491588
## 230  13.9583333 1.3333333  6.502210
## 231  -0.6250000 3.6666667  6.507785
## 232  -1.6666667 2.0000000  6.529267
## 233  -5.8333333 2.8333333  6.533229
## 234   4.5833333 1.0000000  6.560315
## 235   4.5833333 3.0000000  6.563299
## 236   9.7916667 2.1666667  6.567276
## 237  -3.7500000 4.0000000  6.587078
## 238   9.7916667 0.3333333  6.628235
## 239  -7.9166667 3.3333333  6.630803
## 240  -2.7083333 2.1666667  6.654106
## 241  -8.9583333 3.6666667  6.671244
## 242   3.5416667 3.1666667  6.673470
## 243  10.8333333 2.0000000  6.678945
## 244   5.6250000 0.8333333  6.694701
## 245  15.0000000 1.1666667  6.778306
## 246  12.9166667 0.0000000  6.778441
## 247  -1.6666667 3.8333333  6.784418
## 248  -3.7500000 2.3333333  6.812673
## 249   2.5000000 3.3333333  6.817851
## 250  -6.8750000 3.0000000  6.820293
## 251  11.8750000 1.8333333  6.824741
## 252   6.6666667 0.6666667  6.862206
## 253  10.8333333 0.1666667  6.923404
## 254 -10.0000000 4.0000000  6.937598
## 255   1.4583333 3.5000000  6.994324
## 256  -8.9583333 3.5000000  7.000656
## 257  12.9166667 1.6666667  7.002532
## 258  -4.7916667 2.5000000  7.002677
## 259   7.7083333 0.5000000  7.060476
## 260  -2.7083333 4.0000000  7.084841
## 261 -10.0000000 3.8333333  7.089844
## 262  -7.9166667 3.1666667  7.130167
## 263   0.4166667 3.6666667  7.200529
## 264  13.9583333 1.5000000  7.209951
## 265  -5.8333333 2.6666667  7.221637
## 266  11.8750000 0.0000000  7.240382
## 267   8.7500000 0.3333333  7.286999
## 268   0.4166667 1.5000000  7.334256
## 269   1.4583333 1.3333333  7.337808
## 270  -0.6250000 1.6666667  7.364019
## 271   2.5000000 1.1666667  7.374626
## 272 -10.0000000 3.6666667  7.385148
## 273   7.7083333 2.6666667  7.406620
## 274   6.6666667 2.8333333  7.421371
## 275   8.7500000 2.5000000  7.424882
## 276  -1.6666667 1.8333333  7.426698
## 277  -0.6250000 3.8333333  7.433993
## 278   3.5416667 1.0000000  7.444219
## 279  15.0000000 1.3333333  7.444523
## 280  -8.9583333 3.3333333  7.460010
## 281  -6.8750000 2.8333333  7.467006
## 282   5.6250000 3.0000000  7.468939
## 283   9.7916667 2.3333333  7.475916
## 284  -2.7083333 2.0000000  7.521468
## 285   9.7916667 0.1666667  7.539228
## 286   4.5833333 0.8333333  7.545678
## 287   4.5833333 3.1666667  7.548704
## 288  10.8333333 2.1666667  7.559058
## 289  -3.7500000 2.1666667  7.647137
## 290   3.5416667 3.3333333  7.659661
## 291  11.8750000 2.0000000  7.673263
## 292   5.6250000 0.6666667  7.677740
## 293  -1.6666667 4.0000000  7.692235
## 294  -7.9166667 3.0000000  7.736272
## 295   2.5000000 3.5000000  7.800478
## 296  -4.7916667 2.3333333  7.802212
## 297 -10.0000000 3.5000000  7.807292
## 298  10.8333333 0.0000000  7.814675
## 299  12.9166667 1.8333333  7.817171
## 300   6.6666667 0.5000000  7.838860
## 301   1.4583333 3.6666667  7.969574
## 302  -5.8333333 2.5000000  7.984979
## 303  13.9583333 1.6666667  7.989176
## 304  -8.9583333 3.1666667  8.027030
## 305   7.7083333 0.3333333  8.027287
## 306   0.4166667 3.8333333  8.165190
## 307  15.0000000 1.5000000  8.187508
## 308  -6.8750000 2.6666667  8.193587
## 309   8.7500000 0.1666667  8.241149
## 310   0.4166667 1.3333333  8.301388
## 311  -0.6250000 1.5000000  8.313925
## 312   1.4583333 1.1666667  8.318313
## 313 -10.0000000 3.3333333  8.337032
## 314  -1.6666667 1.6666667  8.355791
## 315   2.5000000 1.0000000  8.364519
## 316   7.7083333 2.8333333  8.379466
## 317   8.7500000 2.6666667  8.381953
## 318  -0.6250000 4.0000000  8.385471
## 319   6.6666667 3.0000000  8.406149
## 320   9.7916667 2.5000000  8.413585
## 321  -7.9166667 2.8333333  8.426115
## 322  -2.7083333 1.8333333  8.426548
## 323   3.5416667 0.8333333  8.439528
## 324   5.6250000 3.1666667  8.461726
## 325  10.8333333 2.3333333  8.474036
## 326   9.7916667 0.0000000  8.478520
## 327  -3.7500000 2.0000000  8.525478
## 328   4.5833333 0.6666667  8.542579
## 329   4.5833333 3.3333333  8.545634
## 330  11.8750000 2.1666667  8.562694
## 331  -4.7916667 2.1666667  8.651613
## 332   3.5416667 3.5000000  8.657049
## 333   5.6250000 0.5000000  8.672674
## 334  12.9166667 2.0000000  8.678695
## 335  -8.9583333 3.0000000  8.680643
## 336   2.5000000 3.6666667  8.794926
## 337  -5.8333333 2.3333333  8.803785
## 338  13.9583333 1.8333333  8.820962
## 339   6.6666667 0.3333333  8.828616
## 340 -10.0000000 3.1666667  8.955294
## 341   1.4583333 3.8333333  8.958043
## 342  -6.8750000 2.5000000  8.980670
## 343  15.0000000 1.6666667  8.988245
## 344   7.7083333 0.1666667  9.009064
## 345   0.4166667 4.0000000  9.145049
## 346  -7.9166667 2.6666667  9.180840
## 347   8.7500000 0.0000000  9.212578
## 348  -0.6250000 1.3333333  9.281854
## 349   0.4166667 1.1666667  9.282978
## 350  -1.6666667 1.5000000  9.307069
## 351   1.4583333 1.0000000  9.310430
## 352   8.7500000 2.8333333  9.355442
## 353  -2.7083333 1.6666667  9.358410
## 354   2.5000000 0.8333333  9.363981
## 355   7.7083333 3.0000000  9.365456
## 356   9.7916667 2.6666667  9.371574
## 357   6.6666667 3.1666667  9.401533
## 358  -8.9583333 2.8333333  9.402808
## 359  10.8333333 2.5000000  9.413719
## 360  -3.7500000 1.8333333  9.435451
## 361   3.5416667 0.6666667  9.443186
## 362   5.6250000 3.3333333  9.463375
## 363  11.8750000 2.3333333  9.481528
## 364  -4.7916667 2.0000000  9.537568
## 365   4.5833333 0.5000000  9.547405
## 366   4.5833333 3.5000000  9.550481
## 367  12.9166667 2.1666667  9.574458
## 368 -10.0000000 3.0000000  9.645069
## 369   3.5416667 3.6666667  9.662168
## 370  -5.8333333 2.1666667  9.663967
## 371   5.6250000 0.3333333  9.675832
## 372  13.9583333 2.0000000  9.691784
## 373   2.5000000 3.8333333  9.797596
## 374  -6.8750000 2.3333333  9.813710
## 375   6.6666667 0.1666667  9.827517
## 376  15.0000000 1.8333333  9.832635
## 377   1.4583333 4.0000000  9.955795
## 378  -7.9166667 2.5000000  9.985746
## 379   7.7083333 0.0000000 10.001401
## 380  -8.9583333 2.6666667 10.178945
## 381  -0.6250000 1.1666667 10.262707
## 382  -1.6666667 1.3333333 10.274372
## 383   0.4166667 1.0000000 10.274881
## 384  -2.7083333 1.5000000 10.309794
## 385   1.4583333 0.8333333 10.310810
## 386   8.7500000 3.0000000 10.340712
## 387   9.7916667 2.8333333 10.344239
## 388   7.7083333 3.1666667 10.360838
## 389  -3.7500000 1.6666667 10.368731
## 390   2.5000000 0.6666667 10.370245
## 391  10.8333333 2.6666667 10.371394
## 392 -10.0000000 2.8333333 10.392127
## 393   6.6666667 3.3333333 10.404480
## 394  11.8750000 2.5000000 10.421992
## 395  -4.7916667 1.8333333 10.450784
## 396   3.5416667 0.5000000 10.452788
## 397   5.6250000 3.5000000 10.471342
## 398  12.9166667 2.3333333 10.495695
## 399  -5.8333333 2.0000000 10.555414
## 400   4.5833333 0.3333333 10.557894
## 401   4.5833333 3.6666667 10.560984
## 402  13.9583333 2.1666667 10.592020
## 403   3.5416667 3.8333333 10.672833
## 404  -6.8750000 2.1666667 10.681959
## 405   5.6250000 0.1666667 10.684899
## 406  15.0000000 2.0000000 10.710357
## 407   2.5000000 4.0000000 10.806198
## 408  -7.9166667 2.3333333 10.829649
## 409   6.6666667 0.0000000 10.833033
## 410  -8.9583333 2.5000000 10.997633
## 411 -10.0000000 2.6666667 11.184997
## 412  -0.6250000 1.0000000 11.253106
## 413  -1.6666667 1.1666667 11.253568
## 414   0.4166667 0.8333333 11.274378
## 415  -2.7083333 1.3333333 11.275760
## 416   1.4583333 0.6666667 11.317260
## 417  -3.7500000 1.5000000 11.319554
## 418   9.7916667 3.0000000 11.327800
## 419   8.7500000 3.1666667 11.334693
## 420  10.8333333 2.8333333 11.342504
## 421   7.7083333 3.3333333 11.363145
## 422  11.8750000 2.6666667 11.378723
## 423   2.5000000 0.5000000 11.381508
## 424  -4.7916667 1.6666667 11.384703
## 425   6.6666667 3.5000000 11.412995
## 426  12.9166667 2.5000000 11.436253
## 427   3.5416667 0.3333333 11.466764
## 428  -5.8333333 1.8333333 11.470840
## 429   5.6250000 3.6666667 11.483964
## 430  13.9583333 2.3333333 11.514772
## 431   4.5833333 0.1666667 11.572562
## 432   4.5833333 3.8333333 11.575664
## 433  -6.8750000 2.0000000 11.577499
## 434  15.0000000 2.1666667 11.613857
## 435   3.5416667 4.0000000 11.687606
## 436   5.6250000 0.0000000 11.698346
## 437  -7.9166667 2.1666667 11.704118
## 438  -8.9583333 2.3333333 11.850057
## 439 -10.0000000 2.5000000 12.014612
## 440  -1.6666667 1.0000000 12.241803
## 441  -0.6250000 0.8333333 12.250736
## 442  -2.7083333 1.1666667 12.252859
## 443   0.4166667 0.6666667 12.279612
## 444  -3.7500000 1.3333333 12.283849
## 445   9.7916667 3.1666667 12.319647
## 446  10.8333333 3.0000000 12.323875
## 447   1.4583333 0.5000000 12.328294
## 448  -4.7916667 1.5000000 12.334622
## 449   8.7500000 3.3333333 12.335278
## 450  11.8750000 2.8333333 12.347942
## 451   7.7083333 3.5000000 12.370693
## 452  12.9166667 2.6666667 12.391732
## 453   2.5000000 0.3333333 12.396546
## 454  -5.8333333 1.6666667 12.404937
## 455   6.6666667 3.6666667 12.425724
## 456  13.9583333 2.5000000 12.455038
## 457   3.5416667 0.1666667 12.484048
## 458  -6.8750000 1.8333333 12.494463
## 459   5.6250000 3.8333333 12.500110
## 460  15.0000000 2.3333333 12.537563
## 461   4.5833333 0.0000000 12.590399
## 462   4.5833333 4.0000000 12.593509
## 463  -7.9166667 2.0000000 12.602791
## 464  -8.9583333 2.1666667 12.729440
## 465 -10.0000000 2.3333333 12.873870
## 466  -1.6666667 0.8333333 13.237054
## 467  -2.7083333 1.0000000 13.238626
## 468  -0.6250000 0.6666667 13.253963
## 469  -3.7500000 1.1666667 13.258675
## 470   0.4166667 0.5000000 13.289284
## 471  -4.7916667 1.3333333 13.297115
## 472  10.8333333 3.1666667 13.313236
## 473   9.7916667 3.3333333 13.317929
## 474  11.8750000 3.0000000 13.326923
## 475   8.7500000 3.5000000 13.340982
## 476   1.4583333 0.3333333 13.342870
## 477  -5.8333333 1.5000000 13.353788
## 478  12.9166667 2.8333333 13.358932
## 479   7.7083333 3.6666667 13.382299
## 480  13.9583333 2.6666667 13.409133
## 481   2.5000000 0.1666667 13.414502
## 482  -6.8750000 1.6666667 13.428463
## 483   6.6666667 3.8333333 13.441713
## 484  15.0000000 2.5000000 13.477322
## 485   3.5416667 0.0000000 13.503894
## 486   5.6250000 4.0000000 13.518985
## 487  -7.9166667 1.8333333 13.520842
## 488  -8.9583333 2.0000000 13.630566
## 489 -10.0000000 2.1666667 13.757218
## 490  -2.7083333 0.8333333 14.231261
## 491  -1.6666667 0.6666667 14.237848
## 492  -3.7500000 1.0000000 14.241870
## 493  -0.6250000 0.5000000 14.261607
## 494  -4.7916667 1.1666667 14.269635
## 495   0.4166667 0.3333333 14.302452
## 496  10.8333333 3.3333333 14.308932
## 497  11.8750000 3.1666667 14.313664
## 498  -5.8333333 1.3333333 14.314458
## 499   9.7916667 3.5000000 14.321301
## 500  12.9166667 3.0000000 14.335480
## 501   8.7500000 3.6666667 14.350728
## 502   1.4583333 0.1666667 14.360237
## 503  13.9583333 2.8333333 14.374304
## 504  -6.8750000 1.5000000 14.376179
## 505   7.7083333 3.8333333 14.397107
## 506  15.0000000 2.6666667 14.429997
## 507   2.5000000 0.0000000 14.434759
## 508  -7.9166667 1.6666667 14.454581
## 509   6.6666667 4.0000000 14.460276
## 510  -8.9583333 1.8333333 14.549395
## 511 -10.0000000 2.0000000 14.660302
## 512  -2.7083333 0.6666667 15.229421
## 513  -3.7500000 0.8333333 15.231814
## 514  -1.6666667 0.5000000 15.243095
## 515  -4.7916667 1.0000000 15.250266
## 516  -0.6250000 0.3333333 15.272794
## 517  -5.8333333 1.1666667 15.284720
## 518  11.8750000 3.3333333 15.306664
## 519  10.8333333 3.5000000 15.309725
## 520   0.4166667 0.1666667 15.318423
## 521  12.9166667 3.1666667 15.319589
## 522   9.7916667 3.6666667 15.328763
## 523  -6.8750000 1.3333333 15.335068
## 524  13.9583333 3.0000000 15.348461
## 525   8.7500000 3.8333333 15.363719
## 526   1.4583333 0.0000000 15.379842
## 527  15.0000000 2.8333333 15.393190
## 528  -7.9166667 1.5000000 15.401154
## 529   7.7083333 4.0000000 15.414485
## 530  -8.9583333 1.6666667 15.482776
## 531 -10.0000000 1.8333333 15.579690
## 532  -3.7500000 0.6666667 16.227271
## 533  -2.7083333 0.5000000 16.232086
## 534  -4.7916667 0.8333333 16.237538
## 535  -1.6666667 0.3333333 16.251968
## 536  -5.8333333 1.0000000 16.262857
## 537  -0.6250000 0.1666667 16.286863
## 538  -6.8750000 1.1666667 16.303159
## 539  11.8750000 3.5000000 16.304779
## 540  12.9166667 3.3333333 16.309890
## 541  10.8333333 3.6666667 16.314678
## 542  13.9583333 3.1666667 16.329997
## 543   0.4166667 0.0000000 16.336675
## 544   9.7916667 3.8333333 16.339559
## 545  -7.9166667 1.3333333 16.358333
## 546  15.0000000 3.0000000 16.365044
## 547   8.7500000 4.0000000 16.379354
## 548  -8.9583333 1.5000000 16.428229
## 549 -10.0000000 1.6666667 16.512659
## 550  -3.7500000 0.5000000 17.227287
## 551  -4.7916667 0.6666667 17.230309
## 552  -2.7083333 0.3333333 17.238470
## 553  -5.8333333 0.8333333 17.247529
## 554  -1.6666667 0.1666667 17.263831
## 555  -6.8750000 1.0000000 17.278905
## 556  -0.6250000 0.0000000 17.303308
## 557  12.9166667 3.5000000 17.305321
## 558  11.8750000 3.6666667 17.307125
## 559  13.9583333 3.3333333 17.317657
## 560  10.8333333 3.8333333 17.323067
## 561  -7.9166667 1.1666667 17.324360
## 562  15.0000000 3.1666667 17.344103
## 563   9.7916667 4.0000000 17.353106
## 564  -8.9583333 1.3333333 17.383784
## 565 -10.0000000 1.5000000 17.457033
## 566  -4.7916667 0.5000000 18.227680
## 567  -3.7500000 0.3333333 18.231110
## 568  -5.8333333 0.6666667 18.237678
## 569  -2.7083333 0.1666667 18.247958
## 570  -6.8750000 0.8333333 18.261080
## 571  -1.6666667 0.0000000 18.278188
## 572  -7.9166667 1.0000000 18.297835
## 573  12.9166667 3.6666667 18.305043
## 574  13.9583333 3.5000000 18.310449
## 575  11.8750000 3.8333333 18.313008
## 576  15.0000000 3.3333333 18.329214
## 577  10.8333333 4.0000000 18.334325
## 578  -8.9583333 1.1666667 18.347862
## 579 -10.0000000 1.3333333 18.411055
## 580  -4.7916667 0.3333333 19.228937
## 581  -5.8333333 0.5000000 19.232457
## 582  -3.7500000 0.1666667 19.238145
## 583  -6.8750000 0.6666667 19.248698
## 584  -2.7083333 0.0000000 19.260062
## 585  -7.9166667 0.8333333 19.277628
## 586  13.9583333 3.6666667 19.307583
## 587  12.9166667 3.8333333 19.308391
## 588  -8.9583333 1.0000000 19.319189
## 589  15.0000000 3.5000000 19.319449
## 590  11.8750000 4.0000000 19.321873
## 591 -10.0000000 1.1666667 19.373301
## 592  -5.8333333 0.3333333 20.231184
## 593  -4.7916667 0.1666667 20.233502
## 594  -6.8750000 0.5000000 20.240964
## 595  -3.7500000 0.0000000 20.247912
## 596  -7.9166667 0.6666667 20.262823
## 597  -8.9583333 0.8333333 20.296723
## 598  13.9583333 3.8333333 20.308418
## 599  15.0000000 3.6666667 20.314060
## 600  12.9166667 4.0000000 20.314828
## 601 -10.0000000 1.0000000 20.342604
## 602  -5.8333333 0.1666667 21.233302
## 603  -6.8750000 0.3333333 21.237225
## 604  -4.7916667 0.0000000 21.240905
## 605  -7.9166667 0.5000000 21.252670
## 606  -8.9583333 0.6666667 21.279609
## 607  13.9583333 4.0000000 21.312434
## 608  15.0000000 3.8333333 21.312435
## 609 -10.0000000 0.8333333 21.318001
## 610  -6.8750000 0.1666667 22.236946
## 611  -5.8333333 0.0000000 22.238352
## 612  -7.9166667 0.3333333 22.246546
## 613  -8.9583333 0.5000000 22.267138
## 614 -10.0000000 0.6666667 22.298692
## 615  15.0000000 4.0000000 22.314067
## 616  -6.8750000 0.0000000 23.239678
## 617  -7.9166667 0.1666667 23.243936
## 618  -8.9583333 0.3333333 23.258719
## 619 -10.0000000 0.5000000 23.284009
## 620  -7.9166667 0.0000000 24.244405
## 621  -8.9583333 0.1666667 24.253855
## 622 -10.0000000 0.3333333 24.273388
## 623  -8.9583333 0.0000000 25.252125
## 624 -10.0000000 0.1666667 25.266352
## 625 -10.0000000 0.0000000 26.262495
# Visualisasikan
grid %>% ggplot(aes(b1,b2)) + 
  geom_point(data = filter(grid, rank(dist)<=10), size = 4, color = "red") +
  geom_point(aes(color= -dist))

Yang menarik.. saat kita mencoba memplotkan temuan parameter dari grid ini, kita akan melihat sedikit perbedaan di bandingkan model linear sebelumnya.

ggplot(sim1, aes(x,y))+
  geom_point(size =2, color = "grey30") +
  geom_abline(aes(intercept = b1, slope = b2, color = -dist),
data = filter(grid, rank(dist)<=10))

Terlihat lebih menyempit bukan? Kita bisa memperhalus model di atas dengan terus mengulangi grid search berdasar parameter baru yang kita temukan, hingga kita menemukan satu garis yang benar-benar menggambarkan data.

Namun dari pada berlama2 dengan cara repetitif tersebut, kita bisa langsung tembak dengan menggunakan fungsi optim().

best <- optim(c(0,0),measure_dist,data = sim1)
# Hasil keseluruhan dari fungsi optim
best
## $par
## [1] 4.222248 2.051204
## 
## $value
## [1] 2.128181
## 
## $counts
## function gradient 
##       77       NA 
## 
## $convergence
## [1] 0
## 
## $message
## NULL
# Ambil par nya saja
best$par
## [1] 4.222248 2.051204
# Visual
ggplot(sim1, aes(x,y)) +
  geom_point(size = 2, color = "grey30") +
  geom_abline(intercept = best$par[1], slope = best$par[2])

Par di atas adalah pasangan parameter paling optimal untuk menemukan selisih paling kecil yang bisa didapat dengan terus mengulangi metode grid. Fungsi optim() bisa digunakan untuk seluruh family model (bukan hanya linear)

Terkhusus untuk model linear seperti di atas, kita bisa saja sebenarnya langsung menggunakan satu fungsi yang bisa dengan mudah melakukan fitting model bernama fungsi lm().

sim1_mod <- lm(y~x, data = sim1)
sim1_mod
## 
## Call:
## lm(formula = y ~ x, data = sim1)
## 
## Coefficients:
## (Intercept)            x  
##       4.221        2.052

Perhatikan, bahwa hasil best parameter antara optim() dan lm() sama.

Note, dalam kondisi default, optim() akan selalu mencari nilai par yang menghasilkan value minimal. Jika kita ingin mencari nilai maksimal optimal, ubah argumen control optim()

optim(par, fun, control=list(fnscale = -1)) # ubah fnscale jadi negatif

Catatan tambahan :

Untuk membuat model linear yang lebih robust (tangguh) ada baiknya kita menggunakan ukuran selisih yang berbeda. Selain ukuran root-mean-squared (seperti yang ada dalam fungsi measure_dist), kita bisa gunakan metode ukuran absolut mean, dimana kita menghitung rerata dari nilai absolut selisih y_data dan y_model.

abs_distance <- function(mod, data){
  diff <- data$y - model1(mod, data)
  mean(abs(diff))
}

best2 <- optim(c(0,0), abs_distance, data = sim1)
best2
## $par
## [1] 4.364852 2.048918
## 
## $value
## [1] 1.695664
## 
## $counts
## function gradient 
##      119       NA 
## 
## $convergence
## [1] 0
## 
## $message
## NULL
# visual

ggplot(sim1, aes(x,y)) +
  geom_point(size = 2, color = "grey30") +
  geom_abline(intercept = best2$par[1], slope = best2$par[2])

Coba kita compare distancenya

# measure_dist
measure_dist(c(4.222248,2.051204),sim1)
## [1] 2.128181
# abs_distance
abs_distance(c(4.364852,2.048918), sim1)
## [1] 1.695664

Di sini terlihat bahwa selisih rerata dari metode absolute mean difference lebih kecil dari pada root-mean-squared difference. Kita bisa menggunakan parameter ini untuk model yang lebih baik.

Visualising Model

Seperti yang telah kita pelajari ada beberapa elemen dalam model yang penting untuk dipahami: * Nilai prediksi * Residual (selisih antara y_data dengan y_model, nilai yang tidak terprediksi dalam model) Kita sudah tahu prinsip dan garis besar mencari keduanya, dan bagian ini kita akan lebih menekankan kembali bagaimana kita menggambarkan Kedua elemen di atas sehingga klien lebih mudah memahami kapasitas model kita, serta pola yang muncul dari data

Visualisasi nilai prediksi

Memvisualisasikan model akan lebih mudah dilakukan dengan data_grid(); sebuah fungsi yang mengambil prinsip dari expand.grid(), untuk membuat sebuah data frame baru yang menyusun value berdasarkan nilai unik masing-masing variabel.

Kita kembali pada data ‘sim1’. Pertama, kita membuat data_grid hanya dengan variabel penjelas x

grid2 <- sim1 %>% data_grid(x)
grid2
## # A tibble: 10 x 1
##        x
##    <int>
##  1     1
##  2     2
##  3     3
##  4     4
##  5     5
##  6     6
##  7     7
##  8     8
##  9     9
## 10    10

Kemudian, kita buat kolom baru yang menampung nilai-nilai prediksi sesuai value variabel x. Di sini kita bisa menggunakan fungsi add_predictions().

grid2 <- grid2 %>% 
  add_predictions(sim1_mod)
grid2
## # A tibble: 10 x 2
##        x  pred
##    <int> <dbl>
##  1     1  6.27
##  2     2  8.32
##  3     3 10.4 
##  4     4 12.4 
##  5     5 14.5 
##  6     6 16.5 
##  7     7 18.6 
##  8     8 20.6 
##  9     9 22.7 
## 10    10 24.7

Tahu apa yang unik di sini?

add_predictions() memungkinkan kita untuk menjalankan fungsi sim1_mod dalam objek grid2 yang nota bene tidak lagi memiliki variabel y, sekaligus membentuk kolom baru berisi nilai prediksi. Seandainya, kita menggunakan mutate seperti biasa, kita akan menghasilkan error.

Kemudian visualisasikan

ggplot(sim1,aes(x,y)) +
  geom_point(size = 2, color = "grey30") +
  geom_line(aes(y = pred), data= grid2, color = "red", size = 1)

add_residuals()

Selanjutnya adalah memvisualkan residual atau selisih antar y_data dengan y_model. Untuk bagian ini, kita bisa dengan mudah menggunakan add_residuals() yang fungsinya kurang lebih mirip dengan add_predictions().

sim1 <- sim1 %>% add_residuals(sim1_mod)
sim1
## # A tibble: 30 x 3
##        x     y    resid
##    <int> <dbl>    <dbl>
##  1     1  4.20 -2.07   
##  2     1  7.51  1.24   
##  3     1  2.13 -4.15   
##  4     2  8.99  0.665  
##  5     2 10.2   1.92   
##  6     2 11.3   2.97   
##  7     3  7.36 -3.02   
##  8     3 10.5   0.130  
##  9     3 10.5   0.136  
## 10     4 12.4   0.00763
## # ... with 20 more rows
# Visual

ggplot(sim1, aes(resid)) +
  geom_freqpoly(binwidth = 0.5)

Untuk memberi gambaran yang lebih gamblang, kita bisa membuat plot antara residual dengan var x.

ggplot(sim1, aes(x,resid)) +
  geom_point()+
  geom_ref_line(h = 0)

Plot x-resid di atas menunjukkan kalau residual dari model kita memiliki distribusi yang random, yang artinya model sudah mampu menangkap pattern

Model fitting lain

Modelling adalah tentang membaca pattern data dan berusaha menemukkan hubungan antar variabel, dan kemudian menyusun parameter terbaik yang mampu menggambarkan pergerakan dan hubungan data.

Selain lm() untuk model linier fitting, R juga menyediakan fungsi yang bisa fitting untuk pola pergerakan smooth curve yang disebut ‘loess()’

tryloess <- loess(y~x, data = sim1)
tryloess
## Call:
## loess(formula = y ~ x, data = sim1)
## 
## Number of Observations: 30 
## Equivalent Number of Parameters: 4.19 
## Residual Standard Error: 2.262
# coba kita visualisasikan
grid2_loess <-  add_predictions(grid2, tryloess)
grid2_loess
## # A tibble: 10 x 2
##        x  pred
##    <int> <dbl>
##  1     1  5.34
##  2     2  8.27
##  3     3 10.8 
##  4     4 12.8 
##  5     5 14.6 
##  6     6 16.6 
##  7     7 18.7 
##  8     8 20.8 
##  9     9 22.6 
## 10    10 24.0
ggplot(sim1, aes(x,y)) + geom_point(size = 2, color = "grey30") +
  geom_line(aes(y = pred), data = grid2_loess, color = "red")

Karena dataset sim1 sudah cukup baik digambarkan dengan linier lm(), maka garis smooth curve dari loess() tidak terlalu banyak berbeda.

Selain add_predictions(), R juga memiliki fungsi gather_predictions() yang berfungsi untuk menambah 2 kolom prediksi, juga spread_predictions() yang menambahkan satu kolom untuk setiap model.

# gather_predictions()

grid2 %>% gather_predictions(sim1_mod, tryloess)
## # A tibble: 20 x 3
##    model        x  pred
##    <chr>    <int> <dbl>
##  1 sim1_mod     1  6.27
##  2 sim1_mod     2  8.32
##  3 sim1_mod     3 10.4 
##  4 sim1_mod     4 12.4 
##  5 sim1_mod     5 14.5 
##  6 sim1_mod     6 16.5 
##  7 sim1_mod     7 18.6 
##  8 sim1_mod     8 20.6 
##  9 sim1_mod     9 22.7 
## 10 sim1_mod    10 24.7 
## 11 tryloess     1  5.34
## 12 tryloess     2  8.27
## 13 tryloess     3 10.8 
## 14 tryloess     4 12.8 
## 15 tryloess     5 14.6 
## 16 tryloess     6 16.6 
## 17 tryloess     7 18.7 
## 18 tryloess     8 20.8 
## 19 tryloess     9 22.6 
## 20 tryloess    10 24.0
# spread_predictions()

grid2 %>% spread_predictions(sim1_mod, tryloess)
## # A tibble: 10 x 4
##        x  pred sim1_mod tryloess
##    <int> <dbl>    <dbl>    <dbl>
##  1     1  6.27     6.27     5.34
##  2     2  8.32     8.32     8.27
##  3     3 10.4     10.4     10.8 
##  4     4 12.4     12.4     12.8 
##  5     5 14.5     14.5     14.6 
##  6     6 16.5     16.5     16.6 
##  7     7 18.6     18.6     18.7 
##  8     8 20.6     20.6     20.8 
##  9     9 22.7     22.7     22.6 
## 10    10 24.7     24.7     24.0

Kedua fungsi di atas, berguna jika kita memiliki beberapa model untuk dibandingkan nilai predictor dan residualnya.

Formula dan model famili

Kita sudah belajar dasar tentang model dengan variabel predictor (penjelas) bertipe continuu. Sekarang, bagaimana jika data yang ada adalah categorical?

Prinsipnya, R akan mengonvert seluruh value categorical tersebut menjadi angka. Contoh :

# Untuk kasus ini, kita akan menggunakan dataset 'sim2' yang berada dalam paket modelr
str(sim2)
## tibble [40 x 2] (S3: tbl_df/tbl/data.frame)
##  $ x: chr [1:40] "a" "a" "a" "a" ...
##  $ y: num [1:40] 1.94 1.18 1.24 2.62 1.11 ...
# Kita bisa lihat bahwa sim2 memiliki 2 variabel 'x' dan 'y' di mana 'x' berisi data character

# Coba kita visual

ggplot(sim2, aes(x,y)) + 
  geom_point()

Nah sebelum lebih jauh, perlu kita ingat bahwa model bukanlah alat maha tahu yang akan memberitahu kita soal realita. Kemampuan model tak jauh berbeda dengan kemampuan kita untuk mempresepsikan data itu sendiri, dan matematika yang kita gunakan untuk membangunnya. Oleh karena itu, apa yang kita lakukan setelah ini, hanya sekedar contoh, jangan diambil pusing.

mod2 <- lm(y~x, data = sim2)

grid3 <- sim2 %>% 
  data_grid(x) %>% 
  add_predictions(mod2)

grid3
## # A tibble: 4 x 2
##   x      pred
##   <chr> <dbl>
## 1 a      1.15
## 2 b      8.12
## 3 c      6.13
## 4 d      1.91

Nilai prediksi apa itu? model linear lm() dalam sim2 memprediksi berapa nilai rata-rata y dari masing-masing kategori.

Bagaimana jika sebuah data set memiliki variabel predictor continuu dan kategorikal?

# Untuk contoh kasus ini, kita gunakan sim3
str(sim3)
## tibble [120 x 5] (S3: tbl_df/tbl/data.frame)
##  $ x1 : int [1:120] 1 1 1 1 1 1 1 1 1 1 ...
##  $ x2 : Factor w/ 4 levels "a","b","c","d": 1 1 1 2 2 2 3 3 3 4 ...
##  $ rep: int [1:120] 1 2 3 1 2 3 1 2 3 1 ...
##  $ y  : Named num [1:120] -0.571 1.184 2.237 7.437 8.518 ...
##   ..- attr(*, "names")= chr [1:120] "a" "a" "a" "b" ...
##  $ sd : num [1:120] 2 2 2 2 2 2 2 2 2 2 ...

Kita coba visualkan dulu.

ggplot(sim3) +
  geom_point(aes(x1, y, color = x2))

Ada 2 model yang memungkinkan untuk fit pada data seperti ini :

mod3 <- lm(y~x1 + x2, data = sim3)
mod4 <- lm(y~x1 * x2, data = sim3)

Apa perbedaan diantara keduanya?

untuk mod3, kita akan mengevaluasi efek individual dari masing-masing x1 dan x2. Untuk mod4, kita tidak hanya mengevaluasi efek masing-masing variabel, tapi juga efek dari interaksi x1 dan x2. Secara matematis, mereka akan terlihat seperti ini : y = a_0 + a_1 * x1 + a_2 * x2 + a_12 * x1 * x2

Kita coba visualkan kedua model ini

# Ini modelnya
mod3 <- lm(y~x1 + x2, data = sim3)
mod4 <- lm(y~x1 * x2, data = sim3)

# grid
grid4 <- sim3 %>% 
  data_grid(x1,x2) %>%
  gather_predictions(mod3, mod4)

grid4
## # A tibble: 80 x 4
##    model    x1 x2     pred
##    <chr> <int> <fct> <dbl>
##  1 mod3      1 a      1.67
##  2 mod3      1 b      4.56
##  3 mod3      1 c      6.48
##  4 mod3      1 d      4.03
##  5 mod3      2 a      1.48
##  6 mod3      2 b      4.37
##  7 mod3      2 c      6.28
##  8 mod3      2 d      3.84
##  9 mod3      3 a      1.28
## 10 mod3      3 b      4.17
## # ... with 70 more rows

Setelah itu

ggplot(sim3, aes(x1,y, color = x2)) +
         geom_point() +
         geom_line(data = grid4, aes(y = pred)) +
  facet_wrap(~model)

Mana model yang lebih baik? yang ‘+’ atau yang ’*’?

Untuk mengujinya kita bisa melakukan visualisasi pada residul model

sim3 <- sim3 %>% gather_residuals(mod3, mod4)

# Visual

ggplot(sim3, aes(x1, resid, color = x2)) +
  geom_point()+
  facet_grid(model~x2)

Transformation

Seringkali dalam realita kita tidak akan berhadapan dengan model linear, atau di lain kasus, kita merasa model sederhana kita : y = x_0 + x_1*x1 ternyata tidak terlalu berfungsi. Transformasi model merupakan hal lumrah juga saat kita berhadapan dengan model yang tidak fit memenuhi asumsi OLS.

Untuk melakukan transformasi data, kita bisa memasukkan fungsi transformasinya.

log(y)~sqrt(x1) + x2
y ~ x1 + I(x2/2)
y ~ x1 + I(x2 + 1000)
I(y^2)~x1

Seperti yang terlihat di atas, saat kita mencoba melakokan transformasi data dengan ‘+’,‘-’,’*‘,’/‘,’^’, maka kita harus menggunakan I() agar transformasi kita tidak dianggap sebagai ekspresi dari model itu sendiri.

Famili Model lain

Selain model linier di mana y = a_1 + a_2 * x1 + a_3 * x2 + a_4 * x3 + a_n * x(n-1)… Terdapat banyak model famili lain yang bisa digunakan. Model linier meskipun mudah dipahami, memiliki banyak keterbatasan berupa asumsi-asumsi yang mengekang. Sebagai contoh, adanya asumsi distribusi normal untuk model linier. Untuk menghadapi kondisi data yang beraneka ragam, para ahli statistik telah memiliki beberapa piliham model famili lain.