Intro

Berikut merupakan analisis mengenai prediksi pasien yang dinyatakan memiliki penyakit jantung/tidak pada suatu rumah sakit berdasarkan kategori dari beberapa variabel penunjangnya. Data diperoleh melalui https://www.kaggle.com/ronitf/heart-disease-uci . Dalam melakukan prediksi digunakan metode algoritma regresi logistik dan K-Nearest Neighbor.

1. Setup Library

library(dplyr) # data wrangling
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(gtools) # peluang
## Warning: package 'gtools' was built under R version 4.1.2
library(caret) # confusion matrix
## Warning: package 'caret' was built under R version 4.1.2
## Loading required package: ggplot2
## Loading required package: lattice
library(class) # knn
## Warning: package 'class' was built under R version 4.1.2
library(e1071)
## Warning: package 'e1071' was built under R version 4.1.2
## 
## Attaching package: 'e1071'
## The following object is masked from 'package:gtools':
## 
##     permutations

2. Input Data

heart <- read.csv("heart.csv")

glimpse(heart)
## Rows: 303
## Columns: 14
## $ ï..age   <int> 63, 37, 41, 56, 57, 57, 56, 44, 52, 57, 54, 48, 49, 64, 58, 5~
## $ sex      <int> 1, 1, 0, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 0, 0, 0, 1, 0, 1~
## $ cp       <int> 3, 2, 1, 1, 0, 0, 1, 1, 2, 2, 0, 2, 1, 3, 3, 2, 2, 3, 0, 3, 0~
## $ trestbps <int> 145, 130, 130, 120, 120, 140, 140, 120, 172, 150, 140, 130, 1~
## $ chol     <int> 233, 250, 204, 236, 354, 192, 294, 263, 199, 168, 239, 275, 2~
## $ fbs      <int> 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0~
## $ restecg  <int> 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1~
## $ thalach  <int> 150, 187, 172, 178, 163, 148, 153, 173, 162, 174, 160, 139, 1~
## $ exang    <int> 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0~
## $ oldpeak  <dbl> 2.3, 3.5, 1.4, 0.8, 0.6, 0.4, 1.3, 0.0, 0.5, 1.6, 1.2, 0.2, 0~
## $ slope    <int> 0, 0, 2, 2, 2, 1, 1, 2, 2, 2, 2, 2, 2, 1, 2, 1, 2, 0, 2, 2, 1~
## $ ca       <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0~
## $ thal     <int> 1, 2, 2, 2, 2, 1, 2, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3~
## $ target   <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1~

2.1 Deskripsi Variabel

Berdasarkan data yang telah kita input pada objek heart diperoleh informasi sebagai berikut.

age : umur dalam beberapa tahun

sex : (1 = laki-laki; 0 = perempuan)

cp : tipe nyeri yang paling parah

trestbps : melacak tekanan darah(dalam satuan mm Hg saat masuk ke rumah sakit)

chol : kolestoral (mg / dl)

fbs : (gula darah puasa> 120 mg / dl) (1 = benar; 0 = salah)

restecg : mengembalikan hasil elektrokardiografi

thalach : denyut jantung max yang tercapai

exang : exercise induced angina (1 = ya; 0 = tidak)

oldpeak : ST depresi yang disebabkan oleh olahraga relatif terhadap istirahat

slope : kemiringan segmen ST latihan puncak

ca : jumlah pembuluh darah utama (0-3) diwarnai dengan fluoroskopi

thal : 3 = normal; 6 = cacat tetap; 7 = cacat yang dapat dibalik

target : 1 = sakit atau 0 = tidak sakit

Berikut ini gambaran sedikit pada data yang digunakan.

head(heart)
##   ï..age sex cp trestbps chol fbs restecg thalach exang oldpeak slope ca thal
## 1     63   1  3      145  233   1       0     150     0     2.3     0  0    1
## 2     37   1  2      130  250   0       1     187     0     3.5     0  0    2
## 3     41   0  1      130  204   0       0     172     0     1.4     2  0    2
## 4     56   1  1      120  236   0       1     178     0     0.8     2  0    2
## 5     57   0  0      120  354   0       1     163     1     0.6     2  0    2
## 6     57   1  0      140  192   0       1     148     0     0.4     1  0    1
##   target
## 1      1
## 2      1
## 3      1
## 4      1
## 5      1
## 6      1

3. Data Wrangling

Tahap berikutnya kita akan mengubah tipe data beserta nama dari kolom ï..age

# renaming columns
names(heart) <- c("age","sex", "cp", "trestbps", "chol", "fbs", "restecg", "thalach", "exang", "oldpeak", "slope", "ca", "thal", "target")

heart <- heart %>%
  mutate_at(vars(sex, fbs, exang, target), as.factor)

glimpse(heart)
## Rows: 303
## Columns: 14
## $ age      <int> 63, 37, 41, 56, 57, 57, 56, 44, 52, 57, 54, 48, 49, 64, 58, 5~
## $ sex      <fct> 1, 1, 0, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 0, 0, 0, 1, 0, 1~
## $ cp       <int> 3, 2, 1, 1, 0, 0, 1, 1, 2, 2, 0, 2, 1, 3, 3, 2, 2, 3, 0, 3, 0~
## $ trestbps <int> 145, 130, 130, 120, 120, 140, 140, 120, 172, 150, 140, 130, 1~
## $ chol     <int> 233, 250, 204, 236, 354, 192, 294, 263, 199, 168, 239, 275, 2~
## $ fbs      <fct> 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0~
## $ restecg  <int> 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1~
## $ thalach  <int> 150, 187, 172, 178, 163, 148, 153, 173, 162, 174, 160, 139, 1~
## $ exang    <fct> 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0~
## $ oldpeak  <dbl> 2.3, 3.5, 1.4, 0.8, 0.6, 0.4, 1.3, 0.0, 0.5, 1.6, 1.2, 0.2, 0~
## $ slope    <int> 0, 0, 2, 2, 2, 1, 1, 2, 2, 2, 2, 2, 2, 1, 2, 1, 2, 0, 2, 2, 1~
## $ ca       <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0~
## $ thal     <int> 1, 2, 2, 2, 2, 1, 2, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3~
## $ target   <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1~

Cek missing Value

colSums(is.na(heart))
##      age      sex       cp trestbps     chol      fbs  restecg  thalach 
##        0        0        0        0        0        0        0        0 
##    exang  oldpeak    slope       ca     thal   target 
##        0        0        0        0        0        0

Tidak ada missing value pada data ini

4. Pre-processing

Cek distribusi variabel target

plot(heart$target)

Kelas lebih didominasi oleh pasien yang memiliki penyakit jantung

Cek proporsi target

prop.table(table(heart$target))
## 
##         0         1 
## 0.4554455 0.5445545

Proporsi kedua kelas sudah cukup seimbang, sehingga tidak perlu membutuhkan tahapan lanjut dalam menyeimbangkan proporsi kedua kelas.

**Cross Validation*

Berikutnya dilakukan tahapan splitting data menjadi data_train dan data_test. Hal ini bertujuan data_train digunakan sebagai modeling sedangkan data_test kita gunakan sebagai penguji model yang sudah kita buat jika dihadapkan dengan unseen data. Selain itu hal ini dapat digunakan untuk melihat kemampuan model yang kita buat dalam menghadapi unseen data. Rasio jumlah sampel yang digunakan yaitu 80% sampel data train dan 20% sampel data test

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

intrain <- sample(nrow(heart), nrow(heart)*0.8)
data_train <- heart[intrain,]
data_test <- heart[-intrain,]
# re-check class balance
prop.table(table(data_train$target))
## 
##         0         1 
## 0.4586777 0.5413223

Proporsi kedua kelas masih dalam keadaan seimbang

5. Logistic Regression Modeling

Pemodelan regresi logistik menggunakan fungsi glm(). Adapun tahap dalam melakukan pemodelan ini digunakan kita gunakan semua variabel prediktor

model1 <- glm(target~., 
              family = "binomial", data = data_train)

Model Fitting

Untuk mengetahui variabel apa saja yang dapat meningkatkan peluang variabel target, kita akan mencoba melakukan model fitting menggunakan feature selection dengan metode backward

model2 <- step(model1, direction = "backward", trace = F)


summary(model2)
## 
## Call:
## glm(formula = target ~ sex + cp + trestbps + thalach + oldpeak + 
##     ca + thal, family = "binomial", data = data_train)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.3994  -0.5195   0.1934   0.5620   2.3439  
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  2.345455   2.060778   1.138 0.255062    
## sex1        -1.643499   0.456125  -3.603 0.000314 ***
## cp           1.034075   0.192466   5.373 7.75e-08 ***
## trestbps    -0.025767   0.010944  -2.354 0.018554 *  
## thalach      0.029248   0.009495   3.080 0.002067 ** 
## oldpeak     -0.746310   0.203875  -3.661 0.000252 ***
## ca          -0.587791   0.182139  -3.227 0.001250 ** 
## thal        -0.806133   0.309949  -2.601 0.009299 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 333.83  on 241  degrees of freedom
## Residual deviance: 182.26  on 234  degrees of freedom
## AIC: 198.26
## 
## Number of Fisher Scoring iterations: 5

Hasil dari summary menunjukkan bahwa terdapat 7 variabel prediktor yang berpengaruh positif terhadap variabel target kita yaitu: sex1, cp, trestbps, thalach, oldpeak, ca, dan thal

Selanjutnya, kita akan memprediksi probabilitas dari variabel target kita berdasarkan model2 menggunakan fungsi predict()

data_test$sick_prob <- predict(model2, type = "response", data_test)

Jika kita ingin mengetahui bagaimana distribusi hasil prediksi yang telah kita lakukan dapat divisualisasikan dengan plot density

# visualisasi distribusi hasil prediksi
library(ggplot2)

ggplot(data_test, aes(sick_prob))+
  geom_density(lwd=0.5, fill = 4, alpha = 0.5 )+
  labs(title = "Distribusi Hasil Prediksi Peluang Sakit Jantung",
       x = "Probabilitas")+
  theme_minimal()

Interpretasi dari plot density ini menunjukkan hasil prediksi cenderung menuju nilai 1 yang menandakan bahwa pasien dengan penyakit jantung (1)

Model Evaluation

Untuk mengevaluasi model yang telah kita buat menggunakan fungsi confusionMatrix()

# preparing object
data_test$pred_heart <- ifelse(data_test$sick_prob > 0.5 , yes = "1", no = "0")
data_test$pred_heart <- as.factor(data_test$pred_heart)

# confusion matrix
eval_log <- confusionMatrix(data_test$pred_heart, data_test$target, positive = "1")

eval_log
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  0  1
##          0 21  3
##          1  6 31
##                                           
##                Accuracy : 0.8525          
##                  95% CI : (0.7383, 0.9302)
##     No Information Rate : 0.5574          
##     P-Value [Acc > NIR] : 8.993e-07       
##                                           
##                   Kappa : 0.6975          
##                                           
##  Mcnemar's Test P-Value : 0.505           
##                                           
##             Sensitivity : 0.9118          
##             Specificity : 0.7778          
##          Pos Pred Value : 0.8378          
##          Neg Pred Value : 0.8750          
##              Prevalence : 0.5574          
##          Detection Rate : 0.5082          
##    Detection Prevalence : 0.6066          
##       Balanced Accuracy : 0.8448          
##                                           
##        'Positive' Class : 1               
## 

- Re-call/Sensitivity = dari semua data aktual yang positif, seberapa mampu proporsi model saya menebak benar.

- Specificity = dari semua data aktual yang negatif, seberapa mampu proporsi model saya menebak yang benar.

- Accuracy = seberapa mampu model saya menebak dengan benar target Y.

- Precision = dari semua hasil prediksi, seberapa mampu model saya dapat menebak benar kelas positif.

Model Interpretation

# Odds ratio all coefficients
exp(model2$coefficients) %>% 
  data.frame()
##                      .
## (Intercept) 10.4380249
## sex1         0.1933025
## cp           2.8125035
## trestbps     0.9745617
## thalach      1.0296796
## oldpeak      0.4741126
## ca           0.5555532
## thal         0.4465818

Interpretasi model2:

Berdasarkan hasil tersebut ditunjukkan melalui nilai Odds pasien berjenis kelamin laki-laki (sex1) yaitu 0.19 < 1. Artinya, pasien berjenis kelamin laki-laki lebih kecil (19%) untuk tidak terjangkit penyakit jantung jika dibandingkan dengan Perempuan.

6. K-Nearest Neighbor

Cross Validation

Berikutnya kita siapkan data train dan data test yang baru untuk tahap KNN dengan menggunakan rasio jumlah sampel yang sama seperti sebelumnya

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

index <- sample(nrow(heart), nrow(heart)*0.8)
heart_train_knn <- heart[index,]
heart_test_knn <- heart[-index,]

Data Pre-Processing

Karena tiap prediktor memiliki range yang berbeda maka kita perlu melakukan scaling agar range tiap prediktor tidak berbeda jauh

# Ambil data numeric pada data train dan lakukan scaling
heart_train_x_scaled <- heart_train_knn %>% 
  select_if(is.numeric) %>% 
  scale()

# Ambil data numeric pada data test dan lakukan scaling. Gunakan scale, dan center dari data train

heart_test_x_scaled <- heart_test_knn %>% 
  select_if(is.numeric) %>% 
  scale(center = attr(heart_train_x_scaled, "scaled:center"),
        scale = attr(heart_train_x_scaled, "scaled:scale"))

# Ambil target data train
heart_train_y <- heart_train_knn$target

# Ambil target data test
heart_test_y <- heart_test_knn$target

Cari nilai optimum k

sqrt(nrow(heart_train_knn))
## [1] 15.55635

Nilai k yang akan kita gunakan yaitu 16

KNN

Lakukan tahap KNN

heart_knn <- knn(train = heart_train_x_scaled, 
                    test = heart_test_x_scaled, 
                    cl = heart_train_y, 
                    k = 16)

Evaluation

Evaluasi model menggunakan confusionMatrix()

confusionMatrix(heart_knn, 
                reference = heart_test_y, 
                positive = "1")
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  0  1
##          0 19  2
##          1  8 32
##                                           
##                Accuracy : 0.8361          
##                  95% CI : (0.7191, 0.9185)
##     No Information Rate : 0.5574          
##     P-Value [Acc > NIR] : 3.844e-06       
##                                           
##                   Kappa : 0.66            
##                                           
##  Mcnemar's Test P-Value : 0.1138          
##                                           
##             Sensitivity : 0.9412          
##             Specificity : 0.7037          
##          Pos Pred Value : 0.8000          
##          Neg Pred Value : 0.9048          
##              Prevalence : 0.5574          
##          Detection Rate : 0.5246          
##    Detection Prevalence : 0.6557          
##       Balanced Accuracy : 0.8224          
##                                           
##        'Positive' Class : 1               
## 

Berdasarkan uji coba menggunakan KNN dan regresi logistik, kemampuan regresi logistik dalam memprediksi pasien dengan penyakit jantung lebih baik dibandingkan KNN dikarenakan nilai precision (Pos Pred Value) regresi logistik 84% lebih besar dibandingkan KNN

7. Kesimpulan

Pada kasus ini jika saya diibaratkan sebagai Dokter Penyakit Jantung akan mempertimbangkan treatment yang akan saya lakukan kepada pasien yang terjangkit penyakit jantung/tidak terjangkit penyakit jantung. Dikarenakan treatment di antara dua kategori tersebut berbeda. Apabila treatment yang diberikan kepada pasien yang actual nya tidak sakit jantung akan lebih membahayakan bagi mereka. Oleh sebab itu, saya memutuskan menggunakan metrics precision yang mana saya tidak ingin model prediksi yang saya buat salah memprediksi mana pasien yang benar-benar sakit jantung atau yang tidak sakit jantung.