library(dplyr)
library(data.table)
library(xgboost)
library(Matrix)
library(caret)


##Top categorical and numeric features as indicated by tree models
category_features<-c("L3_S32_F3854", "L3_S32_F3851", "L1_S24_F1523", "L1_S24_F1525", "L2_S26_F3038", "L1_S24_F1530",
                     "L1_S24_F819",  "L2_S27_F3131", "L3_S29_F3475", "L3_S29_F3317")

numeric_features<-c("L3_S32_F3850", "L3_S29_F3412", "L1_S24_F1723", "L3_S33_F3865", "L1_S24_F1846", "L1_S24_F1695",
                    "L3_S33_F3859", "L3_S30_F3519", "L3_S33_F3873", "L3_S29_F3407", "L1_S24_F1632", "L3_S30_F3724",
                    "L3_S30_F3809", "L3_S29_F3382", "L1_S24_F1604", "L3_S30_F3704", "L1_S24_F1844", "L1_S24_F1498",
                    "L3_S38_F3956", "L3_S38_F3960", "L3_S33_F3867", "L3_S33_F3857", "L3_S29_F3351", "L3_S38_F3952",
                    "L1_S24_F867",  "L3_S29_F3339", "L3_S33_F3855", "L0_S0_F20",    "L3_S30_F3744", "L3_S30_F3554",
                    "L3_S34_F3880", "L2_S26_F3040", "L3_S29_F3395", "L3_S30_F3754", "L3_S29_F3479", "L3_S29_F3458",
                    "L2_S27_F3166", "L3_S30_F3804", "L1_S24_F1609", "L3_S30_F3759", "L1_S24_F1842", "L0_S0_F22",   
                    "L3_S36_F3920", "L0_S1_F24",    "L0_S4_F104",   "L3_S29_F3348", "L3_S33_F3861", "L1_S24_F1758",
                    "L0_S12_F350",  "L3_S29_F3336", "L3_S34_F3876", "L3_S30_F3769", "L3_S30_F3774", "L1_S24_F1578",
                    "L0_S1_F28")

##Load in Data
suppressWarnings(dt_cat <- fread("../input/train_categorical.csv",
                                 select = c(category_features),
                                 showProgress = T))

dt_num <- fread("../input/train_numeric.csv",
                select = c(numeric_features, "Response"),
                showProgress = T)

suppressWarnings(dt_cat_test <- fread("../input/test_categorical.csv",
                                      select = c(category_features),
                                      showProgress = T))


dt_num_test <- fread("../input/test_numeric.csv",
                     select = c(numeric_features, "Id"),
                     showProgress = T)


Y<-dt_num$Response
dt_num[ , Response := NULL]

Id  <- dt_num_test$Id
dt_num_test[ , Id := NULL]

row.train=nrow(dt_num)

D.cat<-rbind(dt_cat,dt_cat_test)
D.num<-rbind(dt_num,dt_num_test)

remove(dt_cat,dt_cat_test,dt_num,dt_num_test)
gc(verbose=FALSE)

##Convert columns set NA to 0 and add +2 offset to values
for(col in names(D.num)) set(D.num, j = col, value = D.num[[col]] + 2)
for(col in names(D.num)) set(D.num, which(is.na(D.num[[col]])), col, 0)

#convert categories to numeric integers
for (f in category_features){
  if (class(D.cat[[f]])=='character'){
    levels <- unique(c(D.cat[[f]]))
    D.cat[[f]]<-as.numeric(factor(D.cat[[f]],levels = levels))
  }
}

#Seperate test and training
train.cat=D.cat[1:row.train,]
test.cat=D.cat[(row.train+1):nrow(D.cat),]

remove(D.cat)

train.num=D.num[1:row.train,]
test.num=D.num[(row.train+1):nrow(D.num),]
remove(D.num)

gc(verbose=FALSE)

#combine categorical and numeric features for test and train
train<-cbind(train.num,train.cat)
remove(train.num,train.cat)

test<-cbind(test.num,test.cat)

remove(test.num,test.cat)
gc(verbose=FALSE)

X.train<-  Matrix(as.matrix(train), sparse = T)
X.test<-  Matrix(as.matrix(test), sparse = T)

remove(test,train)
gc(verbose=FALSE)

##create folds or CV
set.seed(7579)
folds <- createFolds(as.factor(Y), k = 6)
valid <- folds$Fold3
model <- c(1:length(Y))[-valid]

##Parameters
param <- list(objective = "binary:logistic",
              eval_metric = "auc",
              eta = 0.02, 
              max_depth = 10, 
              min_child_weight = 0.625, 
              base_score = 0.005)

dmodel <- xgb.DMatrix(X.train[model,], label = Y[model])
dvalid <- xgb.DMatrix(X.train[valid,], label = Y[valid])
remove(X.train)
gc(verbose=FALSE)


##Train Model
set.seed(7579)

m1 <- xgb.train(data = dmodel, param, nrounds = 21,
                watchlist = list(mod = dmodel, val = dvalid), early_stopping_rounds=10)


pred <- predict(m1, dvalid)

remove(dmodel, dvalid)
gc(verbose=FALSE)

#####################################
##Matthews Coefficient to determine threshold

mc <- function(actual, predicted) {
  
  tp <- as.numeric(sum(actual == 1 & predicted == 1))
  tn <- as.numeric(sum(actual == 0 & predicted == 0))
  fp <- as.numeric(sum(actual == 0 & predicted == 1))
  fn <- as.numeric(sum(actual == 1 & predicted == 0))
  
  numer <- (tp * tn) - (fp * fn)
  denom <- ((tp + fp) * (tp + fn) * (tn + fp) * (tn + fn)) ^ 0.5
  
  numer / denom
}

matt <- data.table(thresh = seq(0.990, 0.999, by = 0.0001))

matt$scores <- sapply(matt$thresh, FUN =
                        function(x) mc(Y[valid], (pred > quantile(pred, x)) * 1))

print(matt)


best <- matt$thresh[which(matt$scores == max(matt$scores))]

#######################################
##Create Submission

dtest <- xgb.DMatrix(X.test)
pred  <- predict(m1, dtest)

summary(pred)

sub   <- data.table(Id = Id,
                    Response = (pred > quantile(pred, best)) * 1)

write.csv(sub, "submission file.csv", row.names = F)