## Evaluacija klasifikacije

# Kao i na prethnodnim casovima, pravimo model koji predvidja da li osoba ima
# dijabetes.
library(mlbench)
data("PimaIndiansDiabetes2")
data = PimaIndiansDiabetes2

# Ova baza ima dosta nedostajucih vrednosti
sapply(data, function(x) sum(is.na(x)))

# Radi jednostavnosti ih izbacujemo, mada nije najbolja praksa kada su
# nedostajanja ovako velikog obima.
data = na.omit(data)

# Odvojicemo podatke na trening i test skup, tako da 75% ide na trening.
# Za zadatke klasifikacije je standardno da se podela na trening i test vrsi
# proporcionalnim stratifikovanim uzorkovanjem, gde su stratumi klase ciljne
# promenljive. Ovo garantuje da je odnos klase na trening i test skupu isti.
set.seed(216)
pos_idx = which(data$diabetes == 'pos')
neg_idx = which(data$diabetes == 'neg')
pos_train_idx = sample(pos_idx, floor(0.75*length(pos_idx)))
neg_train_idx = sample(neg_idx, floor(0.75*length(neg_idx)))
train = data[c(pos_train_idx, neg_train_idx),]
test = data[-c(pos_train_idx, neg_train_idx),]

table(train$diabetes) #375:201 ~ 1.87
table(test$diabetes) #125:67 ~ 1.87

# Sada kad znamo kako radi stratifikovana podela, mozemo u buducnosti da je
# vrsimo ugradjenom funkcijom iz paketa caret.
#library(caret)
#train_idx = createDataPartition(data$diabetes, p=0.75, list = FALSE)

# Napravicemo model sa svim prediktora na trening skupu i izvrsicemo njegovu
# evaluaciju na test skupu.
model <- glm(diabetes ~ ., data = train, family = 'binomial')
summary(model)

# Racunamo predvidjene verovatnoce na test skupu:
predicted_probabilities = predict(model, newdata = test, type = 'response')

# Vrsimo klasifikaciju predvidjanja sa standardnim pragom 0.5:
predicted_classes = as.numeric(predicted_probabilities > 0.5)

# Prikazujemo matricu konfuzije:
confusion_matrix = table(test$diabetes, predicted_classes)
confusion_matrix

# paket Metrics ima automatski prikaz matrice konfuzije koji pise jos neke stvari

# Racunamo tacnost modela
acc = (confusion_matrix[1,1] + confusion_matrix[2,2]) / (confusion_matrix[1,1] + confusion_matrix[1,2] + confusion_matrix[2,1] + confusion_matrix[2,2])
acc #solidna vrednost

# Racunamo odziv modela, posto nam je bitno da za osobe koje stvarno imaju dijabetes, to i predvidimo
sens = confusion_matrix[2,2] / (confusion_matrix[2,1] + confusion_matrix[2,2])
sens #los odziv

# Racunamo i f1-score
prec = confusion_matrix[2,2] / (confusion_matrix[1,2] + confusion_matrix[2,2])
prec

f1 = 2*sens*prec/(sens + prec)
f1

# Iscrtavamo i ROC krivu:
library(pROC)
roc_curve = roc(response = test$diabetes, predictor = predicted_probabilities)
plot(roc_curve)
auc(roc_curve)
# Solidno visoka povrsina ispod krive, upucuje da model kvalitetno prepoznaje
# razlike izmedju dijabeticara i nedijabeticara.

# U ovoj situaciji bismo potencijalno hteli bolji odziv, jer zelimo da za 
# dijabeticare model cesce uvidja da imaju dijabetes. Ovo postizemo tako sto
# promenimo prag klasifikacije. Ako zelimo odziv od bar 80%, vidimo sa ROC
# krive da ce specificnost pasti na oko 70%, sto je u redu. Iz podataka za ROC
# krivu izvlacimo najvece C koje nam to garantuje
new_C = roc_curve$thresholds[max(which(roc_curve$sensitivities > 0.8))]
new_C

# Sada racunamo metrike za novi prag. ROC krivu nema potrebe ponovo da pravimo,
# jer ne zavisi od praga/radi sa verovatnocama.

new_predicted_classes = as.numeric(predicted_probabilities > new_C)
new_conf_matrix = table(test$diabetes, new_predicted_classes)
new_conf_matrix

acc = (new_conf_matrix[1,1] + new_conf_matrix[2,2]) / (new_conf_matrix[1,1] + new_conf_matrix[1,2] + new_conf_matrix[2,1] + new_conf_matrix[2,2])
acc #tacnost je opala

sens = new_conf_matrix[2,2] / (new_conf_matrix[2,1] + new_conf_matrix[2,2])
sens #82%

prec = new_conf_matrix[2,2] / (new_conf_matrix[1,2] + new_conf_matrix[2,2])
f1 = 2*sens*prec/(sens + prec)
f1 #f1 malo bolji