1 Objective

Hello everyone, pada kesempatan ini saya akan mencoba membuat model dari data wholesale dimana dalam data tersebut terdapat serangkaian variabel berkaitan dengan Industri FPMG.

Berikut adalah sedikit overview yang bisa kita lihat mengenai apa yang akan kita kerjakan dalam report kali ini

  1. Pemilihan target variabel tergantung dari perspektif case yang ingin diambil
  2. Analisis data dan proses pemilihan variabel prediktor / feature selection
  3. Pre-processing data dari data cleansing hingga cross validation.
  4. Penjelasan mengenai evaluasi dari model yang digunakan, apakah metrik terbaik untuk mengevaluasi model dan mengapa.
  5. Mendokumentasikan analisa tentang cara meningkatkan performa dari modelnya (Misal proses pemilihan k optimum pada model knn). dan/atau perbandingan dari model logistic dan knn.

2 Business Objective

Suatu industri FMCG mempunyai segmentasi klien yaitu Horeca (Hotel, Restaurant, Cafe) & Retail. CRM System mereka telah menyimpan data pembelian dari tiap klien tersebut. Ingin dilakukan klasifikasi segmen klien secara otomatis berdasarkan data pembelian klien tersebut.

Variable target: horeca(1)/retail(2)

3 Data Preparation and Set up

3.1 Library yang kemungkinan akan digunakan

library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.0 --
## v ggplot2 3.3.2     v purrr   0.3.4
## v tibble  3.0.4     v dplyr   1.0.2
## v tidyr   1.1.2     v stringr 1.4.0
## v readr   1.4.0     v forcats 0.5.0
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(lubridate)
## 
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
## 
##     date, intersect, setdiff, union
library(GGally)
## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2
library(MLmetrics)
## 
## Attaching package: 'MLmetrics'
## The following object is masked from 'package:base':
## 
##     Recall
library(lmtest)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(car)
## Loading required package: carData
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
## The following object is masked from 'package:purrr':
## 
##     some
library(plotly)
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
library(caret)
## Loading required package: lattice
## 
## Attaching package: 'caret'
## The following objects are masked from 'package:MLmetrics':
## 
##     MAE, RMSE
## The following object is masked from 'package:purrr':
## 
##     lift
library(class)

3.2 Importing Data

wholesale <- read.csv("data_input/wholesale.csv")

glimpse(wholesale)
## Rows: 440
## Columns: 8
## $ Channel          <int> 2, 2, 2, 1, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 1, 2,...
## $ Region           <int> 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,...
## $ Fresh            <int> 12669, 7057, 6353, 13265, 22615, 9413, 12126, 7579...
## $ Milk             <int> 9656, 9810, 8808, 1196, 5410, 8259, 3199, 4956, 36...
## $ Grocery          <int> 7561, 9568, 7684, 4221, 7198, 5126, 6975, 9426, 61...
## $ Frozen           <int> 214, 1762, 2405, 6404, 3915, 666, 480, 1669, 425, ...
## $ Detergents_Paper <int> 2674, 3293, 3516, 507, 1777, 1795, 3140, 3321, 171...
## $ Delicassen       <int> 1338, 1776, 7844, 1788, 5185, 1451, 545, 2566, 750...

3.2.1 Informasi yang dapat di ambil dari data yang akan kita gunakan

  • Channel: horeca (1), retail (2)
  • Region: lokasi toko; Lisbon (1), Oporto (2), Other regions (3)
  • Fresh: jumlah pembelian produk segar
  • Milk: jumlah pembelian produk susu
  • Grocery: jumlah pembelian produk grocery
  • Frozen: jumlah pembelian produk es/beku
  • Detergents_Paper: jumlah pembelian produk detergent & kertas
  • Delicassen: jumlah pembelian produk delicatessen (berkualitas tinggi)

3.3 Data Wrangling

Untuk data wrangling dikarenakan tujuan kita adalah mengetahui prediski dari variabel channel melalui pemodelan kita maka kita melakukan perubahan tipe data pada kolom variabel ini.

# mengganti nilai di kolom `Channel`
wholesale$Channel <- factor(wholesale$Channel, levels = c(1, 2), labels = c("horeca", "retail"))

head(wholesale)

Selanjutnya kita melihat pada data wholesale masih terdapat kolom Region yang mana dianggap kurang signifikan dalam memberikan pengaruh pada model dikarenakan angka yang relatif sama. maka dari itu kolom variabel ini kita putuskan untuk di eliminasi

wholesale <- wholesale %>%
             select (-Region)

head(wholesale)

3.4 Exploratory Data Analysis

3.4.1 Cek Proporsi Target

prop.table(table(wholesale$Channel))
## 
##    horeca    retail 
## 0.6772727 0.3227273

berdasarkan chunk diatas data yang digunakan dapat di artikan sebagai data balance. inbalance = 90/10 atau nilai perbandingan ekstreme lainnya

Mengapa kita harus melakukan cek proporsi target? …

3.4.2 Cek range nilai antar prediktor (untuk k-NN)

summary(wholesale)
##    Channel        Fresh             Milk          Grocery     
##  horeca:298   Min.   :     3   Min.   :   55   Min.   :    3  
##  retail:142   1st Qu.:  3128   1st Qu.: 1533   1st Qu.: 2153  
##               Median :  8504   Median : 3627   Median : 4756  
##               Mean   : 12000   Mean   : 5796   Mean   : 7951  
##               3rd Qu.: 16934   3rd Qu.: 7190   3rd Qu.:10656  
##               Max.   :112151   Max.   :73498   Max.   :92780  
##      Frozen        Detergents_Paper    Delicassen     
##  Min.   :   25.0   Min.   :    3.0   Min.   :    3.0  
##  1st Qu.:  742.2   1st Qu.:  256.8   1st Qu.:  408.2  
##  Median : 1526.0   Median :  816.5   Median :  965.5  
##  Mean   : 3071.9   Mean   : 2881.5   Mean   : 1524.9  
##  3rd Qu.: 3554.2   3rd Qu.: 3922.0   3rd Qu.: 1820.2  
##  Max.   :60869.0   Max.   :40827.0   Max.   :47943.0

3.5 Cross Validation

Disini kita akan melakukan cross validation dimana data kita akan kita bagi menjadi data trainn dan data test sehingga nantinya kita mampu melakukan validasi pada model kita

RNGkind(sample.kind = "Rounding")
## Warning in RNGkind(sample.kind = "Rounding"): non-uniform 'Rounding' sampler
## used
set.seed(100)

intrain <- sample(nrow(wholesale), nrow(wholesale)* 0.8)
wholesale.train <- wholesale[intrain, ]
wholesale.test <- wholesale[-intrain, ]
# re-check target proportion

prop.table(table(wholesale.train$Channel))
## 
##    horeca    retail 
## 0.6846591 0.3153409

4 Modeling with Logistic Regression

Melakukan pemodelan menggunakan regresi logistik. Pemodelan menggunakan fungsi glm() dalam memodelkan menggunakan regresi logistik. Dalam hal ini seluruh variabel kita anggap mempengaruhi target variabel, dimana variabel Channel menjadi variabel responnya.

model_log <- glm(Channel ~ . , data = wholesale.train, family = "binomial")

4.1 Model Fitting

Berdasarkan pemodelan diatas kita belum mendapatkan variabel mana saja yang memiliki pengaruh paling signifikan dalam pemodelan ini oleh karenanya kita akan melakukan stepwise backward untuk memlih variabel-variabel mana saja yang paling signifikan dalam memberikan input dalam model ini

backward <- step(object = model_log, direction = "backward", trace = F)
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
summary(backward)
## 
## Call:
## glm(formula = Channel ~ Grocery + Frozen + Detergents_Paper, 
##     family = "binomial", data = wholesale.train)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -2.72423  -0.32323  -0.24230   0.04911   3.07108  
## 
## Coefficients:
##                    Estimate Std. Error z value Pr(>|z|)    
## (Intercept)      -3.555e+00  4.338e-01  -8.195 2.50e-16 ***
## Grocery           1.617e-04  5.296e-05   3.054  0.00226 ** 
## Frozen           -1.595e-04  8.050e-05  -1.981  0.04756 *  
## Detergents_Paper  7.489e-04  1.340e-04   5.589 2.28e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 438.81  on 351  degrees of freedom
## Residual deviance: 170.93  on 348  degrees of freedom
## AIC: 178.93
## 
## Number of Fisher Scoring iterations: 7

berdasarkan model yang sudah melalui stepwise backward maka dipilih tiga variabel prediksi dalam case ini yaitu grocery,frozen dan detergent_paper

#Interpretasi dari hasil summary untuk masing2 variabel


#Odds untuk variabel grocery

exp(0.00016174)
## [1] 1.000162
#Odds untuk variabel frozen

exp(-0.00015950)
## [1] 0.9998405
#Odds untuk variabel detergent_paper

exp(0.00074890 )
## [1] 1.000749

4.2 Predict using logistic regression

Dalam tahap ini kita akan melakukan proses menentukan class dari masing-masing hasi peluang yang telah dihasilkan menggunakan fungsi predict()

options(scipen = 123)

predict(backward, wholesale.test, type = "response")
##             1             7            10            18            21 
## 0.40997794155 0.46200433860 0.99241631854 0.05031892770 0.21596549327 
##            23            24            29            31            33 
## 0.07242474024 0.91922643272 0.99899497114 0.41486180758 0.07991103917 
##            35            38            39            41            42 
## 0.05237078509 0.95385567378 0.98727280934 0.03090772272 0.27132285866 
##            46            47            52            56            71 
## 0.99426927901 0.99196671463 0.21245710744 0.19802107845 0.00948802476 
##            84            86            93           100           101 
## 0.03822967250 1.00000000000 0.99999976625 0.01825226789 0.98030488498 
##           103           104           105           111           118 
## 0.66874635026 0.02019283508 0.10961416134 0.02880798884 0.06244232174 
##           123           124           130           137           145 
## 0.02157058834 0.26895385263 0.03632856322 0.21156478411 0.05384732199 
##           149           150           155           158           162 
## 0.01979963512 0.03451144459 0.02819625179 0.04365968772 0.02456508296 
##           165           187           189           193           199 
## 0.54470724141 0.03696798616 0.95332228442 0.04154574038 0.07103548005 
##           208           210           215           223           231 
## 0.71924131566 0.99975647099 0.93074507099 0.01634485015 0.03699677029 
##           236           238           241           247           248 
## 0.08453582813 0.01510301073 0.01831556062 0.03807444204 0.02243064691 
##           251           261           263           265           278 
## 0.03340340480 0.06449233404 0.04230943593 0.96805809358 0.00340723348 
##           284           293           296           299           305 
## 0.00958438774 0.05470484900 0.24567886988 0.78682250510 0.99989664879 
##           306           320           326           332           337 
## 0.66292013705 0.99977762979 0.00004077069 0.99958327383 0.03440231610 
##           338           344           345           363           367 
## 0.01621937585 0.99999417569 0.03035216283 0.04253670995 0.03532828600 
##           372           395           399           401           403 
## 0.01785017471 0.04114725381 0.01384957418 0.02164798585 0.05954776268 
##           404           409           412           415           416 
## 0.08645562106 0.32362320780 0.04337316892 0.04485053824 0.55210879090 
##           418           428           433 
## 0.68664877481 0.00854252731 0.13766100142
log.Risk <- predict(backward, wholesale.test, type = "response")

Menentukan kelas dari hasil peluang yang telah kita dapatkan

# tentukan kelas
log.Label <- ifelse(log.Risk >0.5, "retail", "horeca")

menguubah label ke dalam tipe faktor

log.Label <- as.factor(log.Label)


head(log.Label)
##      1      7     10     18     21     23 
## horeca horeca retail horeca horeca horeca 
## Levels: horeca retail

berdasarkan hasil diatas kita sudah dapat melihat prediksi penentuan class untuk masing-masing row

4.2.1 prediction detail and input prdiction result into our data

summary(backward)
## 
## Call:
## glm(formula = Channel ~ Grocery + Frozen + Detergents_Paper, 
##     family = "binomial", data = wholesale.train)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -2.72423  -0.32323  -0.24230   0.04911   3.07108  
## 
## Coefficients:
##                     Estimate  Std. Error z value            Pr(>|z|)    
## (Intercept)      -3.55540313  0.43384614  -8.195 0.00000000000000025 ***
## Grocery           0.00016174  0.00005296   3.054             0.00226 ** 
## Frozen           -0.00015950  0.00008050  -1.981             0.04756 *  
## Detergents_Paper  0.00074890  0.00013400   5.589 0.00000002283696859 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 438.81  on 351  degrees of freedom
## Residual deviance: 170.93  on 348  degrees of freedom
## AIC: 178.93
## 
## Number of Fisher Scoring iterations: 7

predict(model, newdata, type)

pada type terdapat pilihan:

  • response: menghasilkan peluang
  • link: menghasilkan log of odds

Prediksi probability channel untuk 6 data teratas:

predict(backward, newdata = wholesale.test[1:6,], type = "response")
##          1          7         10         18         21         23 
## 0.40997794 0.46200434 0.99241632 0.05031893 0.21596549 0.07242474

Prediksi probability channel untuk data wholesale.test dan disimpan pada kolom baru bernama pred.Risk.

wholesale.test$pred.Risk <- predict(backward, newdata = wholesale.test, 
                                type = "response")
head(wholesale.test)
# ifelse(kondisi, benar, salah)
wholesale.test$pred.Label <- ifelse(wholesale.test$pred.Risk < 0.5, 
                                yes = "horeca", 
                                no = "retail")
# ubah kelas target (aktual dan prediksi) menjadi factor
wholesale.test$pred.Label <- as.factor(wholesale.test$pred.Label)
# lihat hasil prediksi
wholesale.test %>% 
  select(Channel, pred.Risk, pred.Label) %>% 
  head(6)

5 Modelling using k-NN

wholesale <- read.csv("data_input/wholesale.csv")

glimpse(wholesale)
## Rows: 440
## Columns: 8
## $ Channel          <int> 2, 2, 2, 1, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 1, 2,...
## $ Region           <int> 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,...
## $ Fresh            <int> 12669, 7057, 6353, 13265, 22615, 9413, 12126, 7579...
## $ Milk             <int> 9656, 9810, 8808, 1196, 5410, 8259, 3199, 4956, 36...
## $ Grocery          <int> 7561, 9568, 7684, 4221, 7198, 5126, 6975, 9426, 61...
## $ Frozen           <int> 214, 1762, 2405, 6404, 3915, 666, 480, 1669, 425, ...
## $ Detergents_Paper <int> 2674, 3293, 3516, 507, 1777, 1795, 3140, 3321, 171...
## $ Delicassen       <int> 1338, 1776, 7844, 1788, 5185, 1451, 545, 2566, 750...
RNGkind(sample.kind = "Rounding")
## Warning in RNGkind(sample.kind = "Rounding"): non-uniform 'Rounding' sampler
## used
set.seed(100)

intrain <- sample(nrow(wholesale), nrow(wholesale)* 0.8)
wholesale.train.knn <- wholesale[intrain, ]
wholesale.test.knn <- wholesale[-intrain, ]

Scaling menggunakan: z-score standardization

Pertama-tama, untuk k-NN harus dipisahkan antara data prediktor (x) & data label (y). Kemudian data prediktor akan discaling untuk data train maupun test (berdasarkan mean dan standar deviasi dari data train).

# predictor variables in `train`
wholesale_train_x <- wholesale.train.knn[, -1]

# predictor variables in `test`
wholesale_test_x <- wholesale.test.knn[, -1]

# target variable in `train`
wholesale_train_y <- wholesale.train.knn[, 1]

# target variable in `test`
wholesale_test_y <- wholesale.test.knn[, 1]
# scaling data prediktor (x)

wholesale_train_xs <- scale(data.frame(wholesale_train_x) %>% select_if(is.numeric))



# scale test_x data
wholesale_test_xs <- scale(x = wholesale_test_x, 
                     center = attr(wholesale_train_xs, "scaled:center"), 
                     scale = attr(wholesale_train_xs, "scaled:scale"))
# finding optimum k

sqrt(nrow(wholesale_train_xs))
## [1] 18.76166

Jumlah kelas target: 2 k = 19

library(class)

# k-NN
knn.Label <- knn(train = wholesale_train_xs, 
                 test = wholesale_test_xs,
                 cl = wholesale_train_y,
                 k = 19)

head(knn.Label)
## [1] 2 2 2 1 1 1
## Levels: 1 2

5.1 Model Evaluation Logistic Regresssion

library(caret)

# logistic regression
cm_log <- confusionMatrix(data = log.Label, 
                reference = wholesale.test$Channel)

cm_log
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction horeca retail
##     horeca     56      7
##     retail      1     24
##                                           
##                Accuracy : 0.9091          
##                  95% CI : (0.8287, 0.9599)
##     No Information Rate : 0.6477          
##     P-Value [Acc > NIR] : 0.00000001509   
##                                           
##                   Kappa : 0.7916          
##                                           
##  Mcnemar's Test P-Value : 0.0771          
##                                           
##             Sensitivity : 0.9825          
##             Specificity : 0.7742          
##          Pos Pred Value : 0.8889          
##          Neg Pred Value : 0.9600          
##              Prevalence : 0.6477          
##          Detection Rate : 0.6364          
##    Detection Prevalence : 0.7159          
##       Balanced Accuracy : 0.8783          
##                                           
##        'Positive' Class : horeca          
## 
# k-NN
cm_knn <- confusionMatrix(data = knn.Label  , 
                reference = as.factor(wholesale.test.knn$Channel))

cm_knn
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  1  2
##          1 54  6
##          2  3 25
##                                           
##                Accuracy : 0.8977          
##                  95% CI : (0.8147, 0.9522)
##     No Information Rate : 0.6477          
##     P-Value [Acc > NIR] : 0.00000007518   
##                                           
##                   Kappa : 0.7708          
##                                           
##  Mcnemar's Test P-Value : 0.505           
##                                           
##             Sensitivity : 0.9474          
##             Specificity : 0.8065          
##          Pos Pred Value : 0.9000          
##          Neg Pred Value : 0.8929          
##              Prevalence : 0.6477          
##          Detection Rate : 0.6136          
##    Detection Prevalence : 0.6818          
##       Balanced Accuracy : 0.8769          
##                                           
##        'Positive' Class : 1               
## 

6 Conclusion

berdasarkan hasil summary dari dua model yang ada maka metrics yang digunakan jika ingin memilih sebanyak mungkin kemungkinan yang postif maka kita memilih metric sensitivity namun untuk keakuratan prediksi kita menggunakan metric pos pred value

berdasarkan dua model yang ada maka model yang dipilih adalah model knn dikarenakan ada indikasi perfect seperation di logistic regression