ESS:

# rm(list = ls())
library(ESS)
library(scales)
library(foreach)
library(tidyverse)
library(readxl)
library(tidyverse)
library(glue)

new_data <- genFakeData(
  "rnorm",
  cor_value=0.6,
  nlevel=3,
  n = 100,
  100, 10) %>%
  mutate(Item_ID = row_number(), .before = Loc_RP50)

incon_lv2 <- evalInconsistency(new_data, d_alpha = 1, 2, 15)
incon_lv3 <- evalInconsistency(new_data, d_alpha = 1, 3, 90)
spikePlot(incon_lv2, font_size = 14)
spikePlot(incon_lv3, font_size = 14)

evalInconsistency(new_data, d_alpha = 1, 2, 15)
cut_candi <- new_data$Loc_RP50
cut_info_lv2 <- map_df(cut_candi, ~ evalInconsistency(new_data, d_alpha = 1, 2, .x)$cut_info)
cut_info_lv3 <- map_df(cut_candi, ~ evalInconsistency(new_data, d_alpha = 1, 3, .x)$cut_info)

res = cut_info_lv2 %>%
  rename("Loc_RP50" = "cut_score",
         "L2" = "counts",
         "L2_W" = "weights") %>%
  bind_cols(., cut_info_lv3 %>%
              select(counts, weights) %>%
              rename("L3" = "counts",
                     "L3_W" = "weights"))

alpha_level <- seq(0, 1, by = 0.1)
d1 <- lapply(1:length(alpha_level), function(i) {
# i <- 1
cut_candi <- new_data$Loc_RP50

cut_info_lv2 <- map_df(cut_candi, ~ evalInconsistency(new_data, d_alpha = alpha_level[i], 2, .x)$cut_info)
cut_info_lv3 <- map_df(cut_candi, ~ evalInconsistency(new_data, d_alpha = alpha_level[i], 3, .x)$cut_info)

res = cut_info_lv2 %>%
  rename("Loc_RP50" = "cut_score",
         "L2" = "counts",
         "L2_W" = "weights") %>%
  bind_cols(., cut_info_lv3 %>%
              select(counts, weights) %>%
              rename("L3" = "counts",
                     "L3_W" = "weights"))
WESS = c(which.min(cut_info_lv2$weights),which.min(cut_info_lv3$weights))
CESS = c(which.min(cut_info_lv2$counts),which.min(cut_info_lv3$counts))

WESS_cuts <- res$Loc_RP50[WESS]
CESS_cuts <- res$Loc_RP50[CESS]
c(alpha_level[i], WESS_cuts, CESS_cuts, WESS, CESS)
})

d1 <- do.call('rbind',d1) %>% data.frame() %>% setNames(c("D_alpha","Lv2_W","Lv3_W","Lv2","Lv3","Lv2_cut_weight","Lv3_cut_weight","Lv2_cut_count","Lv3_cut_count")) %>% mutate_all(round, 3)

input <- list(font_size = 12)
d1 %>%
  ggplot() +
  geom_point(aes(x = D_alpha, y = Lv2_W)) +
  geom_line(aes(x = D_alpha, y = Lv2), colour = "red", linetype = "dashed") +

  geom_text(aes(x = D_alpha, y = Lv2_W, label = Lv2_cut_weight),
            hjust=0.5,
            vjust=3,
            size = 5) +

  scale_x_reverse(breaks = seq(0,1, 0.1)) +
  labs(caption = "Numbers below points indicate the cut OOD") +
  theme_bw(base_size = input$font_size) +
  theme(plot.caption=element_text(size=18, hjust=0, face="italic", color="black"))



cut1<- simEstCutScore(new_data, WESS = T, SQRT = F)

res
cut1$res

res$L2_W - cut1$res$L2_W
res$L2
cut1$res$L2

res$L2_W
cut1$res$L2_W


est_data <- list(
res = cut_info_lv2 %>%
  rename("Loc_RP50" = "cut_score",
         "L2" = "counts",
         "L2_W" = "weights") %>%
  bind_cols(., cut_info_lv3 %>%
              select(counts, weights) %>%
              rename("L3" = "counts",
                     "L3_W" = "weights")),

selected_CP = list(
  WESS = c(which.min(cut_info_lv2$weights),which.min(cut_info_lv3$weights)),
  CESS = c(which.min(cut_info_lv2$counts),which.min(cut_info_lv3$counts))
)
)
# which(cut_info_lv2$counts == min(cut_info_lv2$counts))
# which(cut_info_lv2$weights == min(cut_info_lv2$weights))
#
#
# which(cut_info_lv3$counts == min(cut_info_lv3$counts))
# which(cut_info_lv3$weights == min(cut_info_lv3$weights))



mkPlot(.data = est_data, WESS = F, fit = F, font_size = 12)