# What's new?
# Experiment to create clustered profile feature with CLARA object
# Note that I'm using train_sample.csv for training with 100k observations only.
#*****************************************************************
# What's Clara?
# CLARA stands for Clustering Large Applications
# CLARA object means a list representing a clustering of the data into k clusters
# Reference: https://stat.ethz.ch/R-manual/R-devel/library/cluster/html/clara.html

if (!require("pacman")) install.packages("pacman")
p_load(knitr, pryr, caTools, tidyverse, data.table, lubridate, chron, tictoc, viridis,
       cluster, factoextra, DescTools, lightgbm)
set.seed(84)               
options(scipen = 9999, warn = -1, digits= 4)

#*****************************************************************
#Functions + Feature engineering
test_hours <- c("4","5","6","7","8","9","10","11","13","14","15")
most_freq_hours_in_test_data <- c("4","5","9","10","13","14")
least_freq_hours_in_test_data <- c("6","11","15")

#functions
normalize <- function(x) {return((x - min(x)) / (max(x) - min(x)))}

#add features to dataframe
add_features <- function(df) {
    cat(">>>>>Piping it down--->>>>", "\n")
    df <- df %>%
      mutate(day  = day(click_time), 
             hour = hour(click_time),
             hms  = HmsToSec(times(strftime(click_time,"%H:%M:%S"))),
             in_test_hh  = ifelse(hour %in% most_freq_hours_in_test_data, 1,
                           ifelse(hour %in% least_freq_hours_in_test_data, 3, 2))) %>%
      select(-c(click_time)) %>%
      add_count(ip, day, in_test_hh) %>% rename("nip_day_test_hh" = n) %>%
      add_count(os, day, in_test_hh) %>% rename("nos_day_test_hh" = n) %>%
      add_count(app, day, in_test_hh) %>% rename("napp_day_test_hh" = n) %>%
      add_count(channel, day, in_test_hh) %>% rename("nchan_day_test_hh" = n) %>%
      add_count(device, day, in_test_hh) %>% rename("ndev_day_test_hh" = n) %>%
      select(-c(in_test_hh))
      return(df)
  }

#*****************************************************************
# read data and add clustered profile

train_path <- "../input/train_sample.csv"
test_path  <- "../input/test.csv"

vars <- c("ip", "app", "device", "os", "channel", "click_time", "attributed_time", "is_attributed")

total_rows <- 100000
chunk_rows <- 100000
skiprows <- total_rows - chunk_rows

print("preparing data for clustered 💁 profile feature...")
train <- fread(train_path, skip=skiprows, nrows=chunk_rows, col.names = vars, showProgress= FALSE,
               colClasses=list(numeric=1:5)) %>%  select(-c(attributed_time))

test  <- fread(test_path, colClasses=list(numeric=2:6), showProgress= FALSE)
sub <- data.table(click_id = test$click_id, is_attributed = NA)
test$click_id <- NULL
test$is_attributed <- NA

df <- rbind(train, test)

#*****************************************************************
print("creating clustered profile with CLARA...")
#CLARA
profile <- clara(train, 3, samples = 5000, pamLike = TRUE)
# optimal clusters are important. Guesstimated k and sample is just for demonstration.

print("Output from CLARA algorithm")
profile
# fviz_cluster(profile, 
#              palette = viridis(4), 
#              ellipse.type = "t", 
#              geom = "point", pointsize = 1,
#              ggtheme = theme_classic()
#              )

tmp <- cbind(train, cluster_profile = profile$cluster) %>% 
    select(-c(click_time, is_attributed)) %>% distinct()
tmp$cluster_profile <- as.character(tmp$cluster_profile)
df <- left_join(df, tmp, by = c("ip", "app", "device", "os", "channel")) %>% 
    mutate(cluster_profile = if_else(is.na(cluster_profile), "4", cluster_profile))
df$cluster_profile <- as.numeric(df$cluster_profile)
print("Clustered profile added to merged train_test...")
glimpse(df)
rm(profile, tmp)
invisible(gc())
cat("memory in use"); mem_used()
cat("--------------------------------", "\n")

#*****************************************************************

#Split back
print("splitting back to train test...")
train_valid  <- df %>% filter(!is.na(df$is_attributed))
test  <- df %>% filter(is.na(df$is_attributed))
test$is_attributed <- NULL
rm(df)
invisible(gc())

print("Shuffled split for validation...")
sample = sample.split(train_valid$is_attributed, SplitRatio = 0.95)
dtrain = subset(train_valid, sample == TRUE)
dvalid  = subset(train_valid, sample == FALSE)
rm(train_valid, sample)

print("Table of class unbalance")
table(dtrain$is_attributed)
invisible(gc())

cat("train size : ", dim(dtrain), "\n")
cat("valid size : ", dim(dvalid), "\n")

#*****************************************************************
#add features
not_to_normalize <- c("ip", "app", "os", "device", "channel", "day", "hour", "is_attributed", "cluster_profile")

print("Feature engineering...")
tic("Total processing time for train data --->")
dtrain <- add_features(dtrain)
print(object.size(dtrain), units = "auto")
nvars = setdiff(names(dtrain), not_to_normalize)
dtrain[nvars] <- normalize(dtrain[nvars])
print("training data...")
str(dtrain)
toc()

tic("Total processing time for validation data --->")
dvalid <- add_features(dvalid)
print(object.size(dvalid), units = "auto")
nvars = setdiff(names(dvalid), not_to_normalize)
dvalid[nvars] <- normalize(dvalid[nvars])
print("validation data...")
str(dvalid)
toc()

invisible(gc())
cat("memory in use"); mem_used()
cat("--------------------------------", "\n")

categorical_features = c("ip", "app", "os", "device", "channel", "day", "hour")
dtrain = lgb.Dataset(data = as.matrix(dtrain[, colnames(dtrain) != "is_attributed"]), 
                     label = dtrain$is_attributed, categorical_feature = categorical_features)
dvalid = lgb.Dataset(data = as.matrix(dvalid[, colnames(dvalid) != "is_attributed"]), 
                     label = dvalid$is_attributed, categorical_feature = categorical_features)
invisible(gc())
cat("memory in use"); mem_used()
cat("--------------------------------", "\n")

#*****************************************************************
#Modelling

print("Modelling")
params = list(objective = "binary", 
              metric = "auc", 
              learning_rate= 0.01, 
              num_leaves= 7,
              max_depth= 3,
              min_child_samples= 100,
              max_bin= 100,
              subsample= 0.7,
              subsample_freq= 1,
              colsample_bytree= 0.9,
              min_child_weight= 0,
              min_split_gain= 0,
              scale_pos_weight=200
              )

tic("Total time for model training --->")
model <- lgb.train(params, dtrain, valids = list(validation = dvalid), nthread = 4,
                   nrounds = 1000, verbose= 1, early_stopping_rounds = 100, eval_freq = 10)

rm(dtrain, dvalid)
invisible(gc())
toc()
cat("-----------------------------------", "\n")
cat("Validation AUC @ best iter: ", max(unlist(model$record_evals[["validation"]][["auc"]][["eval"]])), "\n")
cat("-----------------------------------", "\n")

#*****************************************************************
#Process test data
print("processing test data...")
tic("Total processing time for test data --->")
test <- add_features(test)
nvars = setdiff(names(test), not_to_normalize)
test[nvars] <- normalize(test[nvars])
print("test data...")
str(test)
print(object.size(test), units = "auto")
test <- as.matrix(test[, colnames(test)])
toc()

cat("memory in use"); mem_used()
cat("--------------------------------", "\n")

print("Predictions")
preds <- predict(model, data = test, n = model$best_iter)
sub$is_attributed = round(as.data.frame(preds),4)
fwrite(sub, "sub_sample_lightgbm_R_clara.csv")
head(sub,15)
cat("-----------------------------------", "\n")

print("Feature importance")
kable(lgb.importance(model, percentage = TRUE))
cat("-----------------------------------", "\n")

print("Preds distribution...")
Freq(sub$is_attributed)
cat("-----------------------------------", "\n")

print("finished...")