Heart Disease Prediction
L. FERUSTA NOVA ADI P.
Kali ini kita akan melakukan prediksi tentang Penyakit Jantung. Dataset yang kita gunakan berasal dari https://archive.ics.uci.edu/ml/datasets/Heart+Disease
Import Library
library(dplyr)
library(tidyr)
library(caret)
library(class)
set.seed(123)Loading Dataset
df <- read.csv("heart.csv", fileEncoding="UTF-8-BOM")str(df)## 'data.frame': 303 obs. of 14 variables:
## $ age : int 63 37 41 56 57 57 56 44 52 57 ...
## $ sex : int 1 1 0 1 0 1 0 1 1 1 ...
## $ cp : int 3 2 1 1 0 0 1 1 2 2 ...
## $ trestbps: int 145 130 130 120 120 140 140 120 172 150 ...
## $ chol : int 233 250 204 236 354 192 294 263 199 168 ...
## $ fbs : int 1 0 0 0 0 0 0 0 1 0 ...
## $ restecg : int 0 1 0 1 1 1 0 1 1 1 ...
## $ thalach : int 150 187 172 178 163 148 153 173 162 174 ...
## $ exang : int 0 0 0 0 1 0 0 0 0 0 ...
## $ oldpeak : num 2.3 3.5 1.4 0.8 0.6 0.4 1.3 0 0.5 1.6 ...
## $ slope : int 0 0 2 2 2 1 1 2 2 2 ...
## $ ca : int 0 0 0 0 0 0 0 0 0 0 ...
## $ thal : int 1 2 2 2 2 1 2 3 3 2 ...
## $ target : int 1 1 1 1 1 1 1 1 1 1 ...Important terms to take note from the data:
age = age in years sex = (1 = male; 0 = female) cp = chest pain type trestbps = resting blood pressure (in mm Hg on admission to the hospital) chol = serum cholestoral in mg/dl fbs = (fasting blood sugar > 120 mg/dl) (1 = true; 0 = false) restecg = resting electrocardiographic results thalach = maximum heart rate achieved exang = exercise induced angina (1 = yes; 0 = no) oldpeak = ST depression induced by exercise relative to rest slope = the slope of the peak exercise ST segment ca = number of major vessels (0-3) colored by flourosopy thal = 3 = normal; 6 = fixed defect; 7 = reversable defect target = 1(heart disease = YES) or 0 (no-heart disease = NO)
Data Pre-Processing
Melihat struktur data di atas, data harus diubah menjadi faktor dari numerik agar sistem dapat menghasilkan output yang sesuai. Kita juga harus memperhitungkan bahwa label 0 dan 1 mungkin tidak mewakili data yang divisualisasikan, oleh karena itu kita harus mengubah levelnya ke data yang lebih dapat dipahami.
heart <- df %<>% mutate(sex = if_else(sex == 1, "MALE", "FEMALE"),
cp = if_else(cp == 1, "ATYPICAL ANGINA",
if_else(cp == 2, "NON-ANGINAL PAIN", "ASYMPTOMATIC")),
fbs = if_else(fbs == 1, ">120", "<=120"),
restecg = if_else(restecg == 0, "NORMAL",
if_else(restecg == 1, "ABNORMALITY", "PROBABLE OR DEFINITE")),
exang = if_else(exang == 1, "YES" ,"NO"),
slope = as.factor(slope),
ca = as.factor(ca),
thal = as.factor(thal),
target = if_else(target == 1, "YES", "NO")) %>%
mutate_if(is.character, 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...
## $ sex <fct> MALE, MALE, FEMALE, MALE, FEMALE, MALE, FEMALE, MALE, MALE...
## $ cp <fct> ASYMPTOMATIC, NON-ANGINAL PAIN, ATYPICAL ANGINA, ATYPICAL ...
## $ trestbps <int> 145, 130, 130, 120, 120, 140, 140, 120, 172, 150, 140, 130...
## $ chol <int> 233, 250, 204, 236, 354, 192, 294, 263, 199, 168, 239, 275...
## $ fbs <fct> >120, <=120, <=120, <=120, <=120, <=120, <=120, <=120, >12...
## $ restecg <fct> NORMAL, ABNORMALITY, NORMAL, ABNORMALITY, ABNORMALITY, ABN...
## $ thalach <int> 150, 187, 172, 178, 163, 148, 153, 173, 162, 174, 160, 139...
## $ exang <fct> NO, NO, NO, NO, YES, NO, NO, NO, NO, NO, NO, NO, NO, YES, ...
## $ 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...
## $ slope <fct> 0, 0, 2, 2, 2, 1, 1, 2, 2, 2, 2, 2, 2, 1, 2, 1, 2, 0, 2, 2...
## $ ca <fct> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2...
## $ thal <fct> 1, 2, 2, 2, 2, 1, 2, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2...
## $ target <fct> YES, YES, YES, YES, YES, YES, YES, YES, YES, YES, YES, YES...Mari kita cek apakah ada NA dalam dataframe
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 0dari tabel di atas dapat terlihat bahwa tidak terdapat NA dalam dataframe
Lalu mari kita lihat proporsi data kita, apakah terdistribusi dengan baik atau tidak
prop.table(table(heart$target))##
## NO YES
## 0.4554455 0.5445545Dapat kita lihat bahwa data terdistribusi cukup baik
Logistic Regression
Cross-validation
Pada part ini kita akan membagi data 80% sebagai data train dan 20% sebagai data test.
set.seed(123)
idx <- sample(nrow(heart), nrow(heart)*0.8)
heart.train <- heart[idx, ]
heart.test <- heart[-idx, ]Modelling
Mari kita analisa data dengan menggunakan model logistic regression (karena target variabelnya merupakan faktor) dengan fungsi glm dengan kolom target sebagai target variabel kita.
model.heart <- glm(target ~ ., heart.train, family="binomial")
summary(model.heart)##
## Call:
## glm(formula = target ~ ., family = "binomial", data = heart.train)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.6615 -0.4079 0.1689 0.4689 2.6053
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.522e+01 2.400e+03 -0.006 0.9949
## age 1.617e-02 2.644e-02 0.612 0.5407
## sexMALE -1.289e+00 5.663e-01 -2.277 0.0228 *
## cpATYPICAL ANGINA 3.890e-01 6.520e-01 0.597 0.5508
## cpNON-ANGINAL PAIN 7.871e-01 5.304e-01 1.484 0.1378
## trestbps -1.551e-02 1.180e-02 -1.314 0.1889
## chol -2.969e-03 4.383e-03 -0.677 0.4982
## fbs>120 6.641e-01 5.869e-01 1.132 0.2578
## restecgNORMAL -7.779e-01 4.326e-01 -1.798 0.0721 .
## restecgPROBABLE OR DEFINITE -1.496e+01 1.620e+03 -0.009 0.9926
## thalach 2.211e-02 1.256e-02 1.760 0.0784 .
## exangYES -1.183e+00 4.668e-01 -2.535 0.0113 *
## oldpeak -3.249e-01 2.534e-01 -1.282 0.1999
## slope1 -9.110e-01 9.246e-01 -0.985 0.3245
## slope2 -3.753e-02 1.047e+00 -0.036 0.9714
## ca1 -1.368e+00 5.445e-01 -2.512 0.0120 *
## ca2 -3.037e+00 7.437e-01 -4.083 4.44e-05 ***
## ca3 -2.197e+00 9.883e-01 -2.223 0.0262 *
## ca4 6.753e-01 1.623e+00 0.416 0.6773
## thal1 1.789e+01 2.400e+03 0.007 0.9941
## thal2 1.755e+01 2.400e+03 0.007 0.9942
## thal3 1.609e+01 2.400e+03 0.007 0.9947
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 333.48 on 241 degrees of freedom
## Residual deviance: 159.90 on 220 degrees of freedom
## AIC: 203.9
##
## Number of Fisher Scoring iterations: 15Setelah itu kita melakukan Backward Stepwise Regression untuk mencari formula yang paling optimal untuk model kita
step(model.heart, direction = "backward")## Start: AIC=203.9
## target ~ age + sex + cp + trestbps + chol + fbs + restecg + thalach +
## exang + oldpeak + slope + ca + thal
##
## Df Deviance AIC
## - cp 2 162.15 202.15
## - age 1 160.27 202.27
## - chol 1 160.35 202.35
## - fbs 1 161.21 203.21
## - slope 2 163.58 203.58
## - oldpeak 1 161.60 203.60
## - trestbps 1 161.65 203.65
## - restecg 2 163.66 203.66
## <none> 159.90 203.90
## - thalach 1 163.14 205.14
## - sex 1 165.25 207.25
## - exang 1 166.44 208.44
## - thal 3 173.12 211.12
## - ca 4 184.99 220.99
##
## Step: AIC=202.15
## target ~ age + sex + trestbps + chol + fbs + restecg + thalach +
## exang + oldpeak + slope + ca + thal
##
## Df Deviance AIC
## - chol 1 162.50 200.50
## - age 1 162.60 200.60
## - oldpeak 1 163.67 201.67
## - fbs 1 163.83 201.83
## - trestbps 1 163.97 201.97
## <none> 162.15 202.15
## - slope 2 166.29 202.29
## - restecg 2 166.95 202.95
## - thalach 1 165.94 203.94
## - sex 1 168.15 206.15
## - exang 1 171.99 209.99
## - thal 3 177.30 211.30
## - ca 4 192.28 224.28
##
## Step: AIC=200.5
## target ~ age + sex + trestbps + fbs + restecg + thalach + exang +
## oldpeak + slope + ca + thal
##
## Df Deviance AIC
## - age 1 162.82 198.82
## - fbs 1 164.10 200.10
## - oldpeak 1 164.14 200.14
## - trestbps 1 164.26 200.26
## <none> 162.50 200.50
## - slope 2 166.88 200.88
## - thalach 1 166.04 202.04
## - restecg 2 168.06 202.06
## - sex 1 168.22 204.22
## - exang 1 172.18 208.18
## - thal 3 178.35 210.35
## - ca 4 192.50 222.50
##
## Step: AIC=198.82
## target ~ sex + trestbps + fbs + restecg + thalach + exang + oldpeak +
## slope + ca + thal
##
## Df Deviance AIC
## - trestbps 1 164.31 198.31
## - fbs 1 164.39 198.39
## - oldpeak 1 164.50 198.50
## <none> 162.82 198.82
## - slope 2 167.09 199.09
## - thalach 1 166.05 200.05
## - restecg 2 168.31 200.31
## - sex 1 168.93 202.93
## - exang 1 172.81 206.81
## - thal 3 178.91 208.91
## - ca 4 194.08 222.08
##
## Step: AIC=198.31
## target ~ sex + fbs + restecg + thalach + exang + oldpeak + slope +
## ca + thal
##
## Df Deviance AIC
## - fbs 1 165.45 197.45
## - oldpeak 1 166.12 198.12
## <none> 164.31 198.31
## - slope 2 168.47 198.47
## - thalach 1 167.25 199.25
## - restecg 2 171.30 201.30
## - sex 1 169.84 201.84
## - exang 1 173.84 205.84
## - thal 3 181.64 209.64
## - ca 4 195.04 221.04
##
## Step: AIC=197.45
## target ~ sex + restecg + thalach + exang + oldpeak + slope +
## ca + thal
##
## Df Deviance AIC
## - slope 2 169.39 197.39
## <none> 165.45 197.45
## - oldpeak 1 167.57 197.57
## - thalach 1 168.84 198.84
## - restecg 2 172.42 200.42
## - sex 1 170.84 200.84
## - exang 1 175.00 205.00
## - thal 3 182.55 208.55
## - ca 4 195.05 219.05
##
## Step: AIC=197.39
## target ~ sex + restecg + thalach + exang + oldpeak + ca + thal
##
## Df Deviance AIC
## <none> 169.39 197.39
## - sex 1 174.43 200.43
## - oldpeak 1 174.44 200.44
## - restecg 2 176.90 200.90
## - thalach 1 175.50 201.50
## - exang 1 179.71 205.71
## - thal 3 189.30 211.30
## - ca 4 198.50 218.50##
## Call: glm(formula = target ~ sex + restecg + thalach + exang + oldpeak +
## ca + thal, family = "binomial", data = heart.train)
##
## Coefficients:
## (Intercept) sexMALE
## -17.7209 -1.0933
## restecgNORMAL restecgPROBABLE OR DEFINITE
## -1.0101 -15.9557
## thalach exangYES
## 0.0261 -1.3795
## oldpeak ca1
## -0.4442 -0.9912
## ca2 ca3
## -2.8911 -1.7628
## ca4 thal1
## 0.6167 17.8994
## thal2 thal3
## 17.6510 15.9149
##
## Degrees of Freedom: 241 Total (i.e. Null); 228 Residual
## Null Deviance: 333.5
## Residual Deviance: 169.4 AIC: 197.4setelah kita menemukan formula yang paling optimal, maka kita akan menyimpannya ke obyek model.backward
model.backward <- glm(formula = target ~ age + sex + cp + trestbps + chol + thalach +
exang + slope + ca + thal, family = "binomial", data = heart.train)
summary(model.backward)##
## Call:
## glm(formula = target ~ age + sex + cp + trestbps + chol + thalach +
## exang + slope + ca + thal, family = "binomial", data = heart.train)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.6556 -0.4409 0.1656 0.4831 2.9197
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -13.172417 882.748163 -0.015 0.98809
## age 0.017502 0.025824 0.678 0.49793
## sexMALE -1.405849 0.532554 -2.640 0.00829 **
## cpATYPICAL ANGINA 0.556275 0.640231 0.869 0.38492
## cpNON-ANGINAL PAIN 0.954302 0.506173 1.885 0.05939 .
## trestbps -0.017546 0.010854 -1.617 0.10598
## chol -0.004529 0.004128 -1.097 0.27261
## thalach 0.024718 0.012104 2.042 0.04113 *
## exangYES -1.174817 0.456570 -2.573 0.01008 *
## slope1 -0.409281 0.778331 -0.526 0.59900
## slope2 0.750759 0.807965 0.929 0.35279
## ca1 -1.421165 0.527395 -2.695 0.00705 **
## ca2 -3.061487 0.684940 -4.470 7.83e-06 ***
## ca3 -2.248931 0.942719 -2.386 0.01705 *
## ca4 1.211613 1.645783 0.736 0.46161
## thal1 14.859017 882.743821 0.017 0.98657
## thal2 14.448510 882.743551 0.016 0.98694
## thal3 13.056795 882.743533 0.015 0.98820
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 333.48 on 241 degrees of freedom
## Residual deviance: 167.69 on 224 degrees of freedom
## AIC: 203.69
##
## Number of Fisher Scoring iterations: 13Prediction
Kita akan menggunakan model yang kita train untuk memprediksi data test kita
pred <- predict(model.backward, heart.test, type = "response")
pred_label <- as.factor(ifelse(pred > 0.5,"YES", "NO"))Model Evaluation
Setelah itu kita akan memprediksi model kita ke data test
library(caret)
confusionMatrix(pred_label, heart.test$target, positive = "YES")## Confusion Matrix and Statistics
##
## Reference
## Prediction NO YES
## NO 22 3
## YES 6 30
##
## Accuracy : 0.8525
## 95% CI : (0.7383, 0.9302)
## No Information Rate : 0.541
## P-Value [Acc > NIR] : 2.6e-07
##
## Kappa : 0.7005
##
## Mcnemar's Test P-Value : 0.505
##
## Sensitivity : 0.9091
## Specificity : 0.7857
## Pos Pred Value : 0.8333
## Neg Pred Value : 0.8800
## Prevalence : 0.5410
## Detection Rate : 0.4918
## Detection Prevalence : 0.5902
## Balanced Accuracy : 0.8474
##
## 'Positive' Class : YES
## Dari confusion matrix di atas bisa kita lihat bahwa model Logistic Regression cukuo bagus dalam memprediksi kelas positif dan negatifnya, akurasi untuk kelas YES dan No tidak berbeda jauh yaitu hanya selisih 0.11%.
K-Nearest Neighbour
Cross Validation
set.seed(100)
idx <- sample(nrow(heart), nrow(heart)*0.8)
train <- heart[idx,]
test <- heart[-idx,]
train_label <- train$target
test_label <- test$targetData Scaling
langkah penting dalam penggunaan KNN adalah Data scaling, setelah Data scaling dilakukan selanjutnya kita akan k-value dengan fungsi sqrt
train_knn <- train %>%
select_if(is.numeric) %>%
scale()
test_knn <-test %>%
select_if(is.numeric) %>%
scale(center = attr(train_knn, "scaled:center"),
scale = attr(train_knn, "scaled:scale"))
sqrt(nrow(train_knn))## [1] 15.55635diketahui bahwa nilai K nya adalah 15.56 dan kita akan bulatkan menjadi 15
Modelling Prediction
library(class)
model_knn <- knn(train = train_knn,
test = test_knn,
cl = train_label,
k = 15)
summary(model_knn)## NO YES
## 20 41Model Evaluation
confusionMatrix(model_knn, test$target, positive = "YES")## Confusion Matrix and Statistics
##
## Reference
## Prediction NO YES
## NO 15 5
## YES 12 29
##
## Accuracy : 0.7213
## 95% CI : (0.5917, 0.8285)
## No Information Rate : 0.5574
## P-Value [Acc > NIR] : 0.00633
##
## Kappa : 0.4197
##
## Mcnemar's Test P-Value : 0.14561
##
## Sensitivity : 0.8529
## Specificity : 0.5556
## Pos Pred Value : 0.7073
## Neg Pred Value : 0.7500
## Prevalence : 0.5574
## Detection Rate : 0.4754
## Detection Prevalence : 0.6721
## Balanced Accuracy : 0.7042
##
## 'Positive' Class : YES
## Dari confusion matrix di atas bisa kita lihat bahwa KNN lebih bagus dalam memprediksi kelas positif saja, aku rasi untuk kelas YES bisa mencapai 90% sedangkan untuk kelas NO sangat rendah yaitu hanya 54.84%.
Conclusion
Dari 2 model evaluation di atas diketahui bahwa model Logistic Regression memiliki angka accuracy yang lebih bagus yaitu di angka 88.52% sedangkan pada model KNN nilai accuracy-nya hanya berada pada angka 72.13%. Metric sensitivity dan specifity juga berbanding lurus dengan accuracy, model KNN memiliki sensitivity dan specifity yang lebih rendah dari model Logistic regression. Menurut saya, karena jumlah dataset yang digunakan terlalu kecil maka hasil di sini belum bisa menjadi patokan dalam menentukan Model mana yang lebih bagus. Untuk dapat mengimprove hasil dari tiap model hal yang perlu dilakukan adalah memperbesar jumlah dataset yang ditrain.
Namun jika merujuk pada hasil dari dataset di atas, pada model regresi kita, variabel atau prediktor yang paling berpengaruh adalah ca yang merupakan jumlah dari major vessels, cp yang merupakan chestpain atau nyeri dada dan pasien berjenis kelamin laki-laki atau variabel sex dengan label ’MALE`.