file.name<-'http://archive.ics.uci.edu/ml/machine-learning-databases/statlog/german/german.data-numeric'
 
credit.approval<-read.table(file.name)
 
names(credit.approval)[25]<-'target'
 
 
######Podział na zbiory
 
 
dziel<-function(dane, frakcje)
{
  frakcja_train<-frakcje[1]
  frakcja_valid<-frakcje[2]
  frakcja_test<-frakcje[3]
 
  losuj<-runif(dim(dane)[1])
 
  assign<-ifelse(losuj<frakcja_train, 'train', ifelse(losuj<frakcja_train+frakcja_valid, 'valid', 'test'))
  assign<-factor(assign, levels=c('train', 'valid', 'test'))
  return(assign)
 
}
 
 
dziel2<-function(dane, frakcje)
{
 
  losuj<-runif(dim(dane)[1])
 
 
  frakcja_train<-quantile(losuj, frakcje[1])
  frakcja_valid<-quantile(losuj, frakcje[2]+frakcje[1])
#   frakcja_test<-quantile(losuj, frakcje[3])
 
 
 
 
 
  assign<-ifelse(losuj<frakcja_train, 'train', ifelse(losuj<frakcja_valid, 'valid', 'test'))
  assign<-factor(assign, levels=c('train', 'valid', 'test'))
  return(assign)
 
}
 
 
 
 
 
 
table(credit.approval$target)
 
credit.approval$target<-credit.approval$target-1
 
przydzial<-dziel2(credit.approval, c(.6,.2,.2))
 
dane$train<-credit.approval[which(przydzial=='train'),]
dane$valid<-credit.approval[which(przydzial=='valid'),]
dane$test<-credit.approval[which(przydzial=='test'),]
 
 
dane<-split(credit.approval, przydzial)
 
 
ols2<-glm(target~., data=dane$train)
logit2<-glm(target~., data=dane$train, family=binomial("logit"))
probit2<-glm(target~., data=dane$train, family=binomial("probit"))
 
ols2_step<-step(ols2)
logit2_step<-step(logit2)
probit2_step<-step(probit2)
 
 
 
predict(ols2, newdata=dane$valid, type='response')
predict(logit2, newdata=dane$valid, type='response')
predict(probit2, newdata=dane$valid, type='response')
 
 
koszt2 <- function(prog, model, data) {
decyzja <- ifelse(predict(model, newdata=data, type='response')>prog,1,0) 
decyzja<-factor(decyzja, levels=c(0,1))
pt <- prop.table(table(decyzja,data$target)) 
pt[1,2] + pt[2,1]
}
 
 
 
 
 
 
 
optymalne_odciecie<-c(optimize(koszt2, c(0,1), model=ols2, data=dane$train)$minimum,
optimize(koszt2, c(0,1), model=ols2_step, data=dane$train)$minimum,
optimize(koszt2, c(0,1), model=logit2, data=dane$train)$minimum,
optimize(koszt2, c(0,1), model=logit2_step, data=dane$train)$minimum,
optimize(koszt2, c(0,1), model=probit2, data=dane$train)$minimum,
optimize(koszt2, c(0,1), model=probit2_step, data=dane$train)$minimum)
 
 
 
koszt_train<-c(
  koszt2(optymalne_odciecie[1], ols2, data=dane$train),
  koszt2(optymalne_odciecie[2], ols2_step, data=dane$train),
  koszt2(optymalne_odciecie[3], logit2, data=dane$train),
  koszt2(optymalne_odciecie[4], logit2_step, data=dane$train),
  koszt2(optymalne_odciecie[5], probit2, data=dane$train),
  koszt2(optymalne_odciecie[6], probit2_step, data=dane$train))
 
 
koszt_valid<-c(
  koszt2(optymalne_odciecie[1], ols2, data=dane$valid),
  koszt2(optymalne_odciecie[2], ols2_step, data=dane$valid),
  koszt2(optymalne_odciecie[3], logit2, data=dane$valid),
  koszt2(optymalne_odciecie[4], logit2_step, data=dane$valid),
  koszt2(optymalne_odciecie[5], probit2, data=dane$valid),
  koszt2(optymalne_odciecie[6], probit2_step, data=dane$valid))
 
 
credit.approval
 
 
 
podsumowanie<-data.frame(model=c("OLS","OLS Step", "Logit", "Logit Step", "Probit", "Probit Step" ),cutoff=optymalne_odciecie, KosztTrain=koszt_train, KosztValid=koszt_valid)
 
 
 
set.seed(38166)
zmienne_od_czapy<-replicate(200,rnorm(dim(dane$train)[1]))
dane_od_czapy<-cbind(dane$train, zmienne_od_czapy)
 
klej4<-function(x)
{paste('V', x)}
 
names(dane_od_czapy)[26:45]<-sapply(26:45, klej4)
 
 
ols_od_czapy<-glm(target~., data=dane_od_czapy)
opt_od_czapy<-optimize(koszt2, c(0,1), model=ols_od_czapy, data=dane_od_czapy)
 
 
dane_od_czapy_valid<-cbind(dane$valid, zmienne_od_czapy)
names(dane_od_czapy_valid)[26:45]<-sapply(26:45, klej4)
 
Z<-koszt2(opt_od_czapy$minimum, ols_od_czapy, data=dane_od_czapy_valid)
 
 
podsumowanie<-rbind(podsumowanie, data.frame(model="od czapy", cutoff=opt_od_czapy$minimum, KosztTrain=opt_od_czapy$objective, KosztValid=Z))
podsumowanie
 
 
 
 
##########KRZYWE LIFT, GAIN, ROC
 
library('ROCR')
 
 
prognoza<-predict(probit2_step, newdata=dane$valid, type="response")
klasa<-ifelse(prognoza>podsumowanie$cutoff[6],1,0)
table(prognoza=klasa, realizacja=dane$valid$target)
 
 
pred <- prediction(prognoza, dane$valid$target)
pred_wizard<- prediction(dane$valid$target,dane$valid$target)
pred_random<-prediction(runif(length(dane$valid$target),0,1), dane$valid$target)
 
 
 
 
 
 
 
 
 
 
#########TRAIN SET
prognoza_train<-predict(probit2_step, newdata=dane$train, type="response")
klasa_train<-ifelse(prognoza_train>podsumowanie$cutoff[6],1,0)
table(klasa_train, dane$train$target)
 
pred_train <- prediction(prognoza_train, dane$train$target)
 
 
 
 
###############GAIN
 
 
 
krzywa_gain <- performance(pred,"tpr","rpp")
krzywa_gain_wizard <- performance(pred_wizard,"tpr","rpp")
krzywa_gain_random <- performance(pred_random,"tpr","rpp")
 
krzywa_gain_train <- performance(pred_train,"tpr","rpp")
 
 
plot(krzywa_gain_train, col="green", lty=2)
plot(krzywa_gain, col="green", add=TRUE)
plot(krzywa_gain_wizard, add=TRUE, col="blue")
plot(krzywa_gain_random, add=TRUE, col="red")
 
 
 
 
 
klasyfikacja <- function(prog, model, data) {
decyzja <- ifelse(predict(model, newdata=data, type='response')>prog,1,0) 
decyzja<-factor(decyzja, levels=c(0,1))
pt <- (table(prognoza=decyzja,realizacja=data$target))
 
return(pt)
}
 
rpp <- function(prog, model, data) {
decyzja <- ifelse(predict(model, newdata=data, type='response')>prog,1,0) 
decyzja<-factor(decyzja, levels=c(0,1))
pt <- prop.table(table(decyzja,data$target)) 
pt[2,1] + pt[2,2]
}
 
 
tpr <- function(prog, model, data) {
decyzja <- ifelse(predict(model, newdata=data, type='response')>prog,1,0) 
decyzja<-factor(decyzja, levels=c(0,1))
pt <- prop.table(table(decyzja,data$target)) 
pt[2,2] / (pt[2,2]+pt[1,2])
}
 
 
 
 
 
 
 
 
##############LIFT
 
krzywa_lift <- performance(pred,"lift","rpp")
krzywa_lift_wizard <- performance(pred_wizard,"lift","rpp")
krzywa_lift_random <- performance(pred_random,"lift","rpp")
 
krzywa_lift_train <- performance(pred_train,"lift","rpp")
 
 
 
plot(krzywa_lift, col="green", lwd=3)
plot(krzywa_lift_wizard, add=TRUE, col="blue")
plot(krzywa_lift_random, add=TRUE, col="red")
plot(krzywa_lift_train, add=TRUE, col="orange", lty=3, lwd=1.5)
grid()
 
 
 
lift <- function(prog, model, data) {
decyzja <- ifelse(predict(model, newdata=data, type='response')>prog,1,0) 
decyzja<-factor(decyzja, levels=c(0,1))
pt <- prop.table(table(decyzja,data$target)) 
(pt[2,2]/(pt[2,2]+pt[1,2])) / (pt[2,2]+pt[2,1])
}
 
 
#############ROC
krzywa_roc <- performance(pred,"tpr","fpr")
krzywa_roc_wizard <- performance(pred_wizard,"tpr","fpr")
krzywa_roc_random <- performance(pred_random,"tpr","fpr")
 
krzywa_roc_train <- performance(pred_train,"tpr","fpr")
 
 
plot(krzywa_roc, col="green")
plot(krzywa_roc_train, col="green", add=TRUE, lty=2)
plot(krzywa_roc_wizard, add=TRUE, col="blue")
plot(krzywa_roc_random, add=TRUE, col="red")
grid()
 
 
fpr <- function(prog, model, data) {
decyzja <- ifelse(predict(model, newdata=data, type='response')>prog,1,0) 
decyzja<-factor(decyzja, levels=c(0,1))
pt <- prop.table(table(decyzja,data$target)) 
pt[2,1] / (pt[2,1]+pt[1,1])
}
 
##############PORÓWNYWANIE MODELI
 
 
 
prognoza_probit2<-predict(probit2, newdata=dane$valid, type="response")
klasa_probit2<-ifelse(prognoza>podsumowanie$cutoff[5],1,0)
pred_probit2 <- prediction(klasa_probit2, dane$valid$target)
lift_probit2<-performance(pred_probit2,"lift","rpp")
 
prognoza_ols2<-predict(ols2, newdata=dane$valid, type="response")
klasa_ols2<-ifelse(prognoza>podsumowanie$cutoff[1],1,0)
pred_ols2 <- prediction(klasa_ols2, dane$valid$target)
lift_ols2<-performance(pred_ols2,"lift","rpp")
 
 
plot(krzywa_lift, col="green")
plot(lift_ols2, col="blue", add=TRUE)
plot(lift_probit2, col="mediumpurple", add=TRUE)
 
 
 
 
krzywa_gain <- performance(pred,"tpr","rpp")
krzywa_probit2 <- performance(pred_probit2,"tpr","rpp")
krzywa_ols2 <- performance(pred_ols2,"tpr","rpp")
 
 
plot(krzywa_gain, col="green")
plot(krzywa_probit2, col="blue", add=TRUE)
plot(krzywa_ols2, col="mediumpurple", add=TRUE)
grid()
 
 
 
krzywa_roc <- performance(pred,"tpr","fpr")
krzywa_roc_probit2 <- performance(pred_probit2,"tpr","fpr")
krzywa_roc_ols2 <- performance(pred_ols2,"tpr","fpr")
 
plot(krzywa_roc, col="green")
plot(krzywa_roc_probit2, col="blue", add=TRUE)
plot(krzywa_roc_ols2, col="mediumpurple", add=TRUE)
grid()