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$targetCari 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.