library(fpca)
library(tidyverse)

# Set working directory.
#setwd("~/.../ACB/Code")

df_pir_acb <- read.csv("df_pir_acb.csv")
df_pir_acb <- df_pir_acb %>%
  select(-Season)

# Candidate models for fitting:
M.set <- c(5, 10, 15, 20)
r.set <- c(3, 4, 5)

#Fit candidate models by fpca.mle:
s0 <- Sys.time()
result_pir_acb <- fpca.mle(as.data.frame(df_pir_acb), M.set, r.set)
s1 <- Sys.time() - s0
s1
#Time difference of 15.05944 mins
#save(result_pir_acb, file = "result_pir_acb.RData")
#load("result_pir_acb.RData")

#Model selection result, selected with the smallest CV score among all converged models:
M <- result_pir_acb$selected_model[1]
r <- result_pir_acb$selected_model[2]

# Look at the CV scores and convergence for each model:
result_pir_acb$cv_scores
result_pir_acb$converge
# Rescaled grid:
grids.new <- result_pir_acb$grid
# Estimated eigenvalues:
evalest <- result_pir_acb$eigenvalues
evalest

# Estimated error variance:
sig2est <- result_pir_acb$error_var

# Estimated eigenfunctions:
eigenfest <- result_pir_acb$eigenfunctions

va <- evalest[1] * eigenfest[1,] * eigenfest[1,] + 
  evalest[2] * eigenfest[2,] * eigenfest[2,] + 
  evalest[3] * eigenfest[3,] * eigenfest[3,] + 
  evalest[4] * eigenfest[4,] * eigenfest[4,] + 
  evalest[5] * eigenfest[5,] * eigenfest[5,] 
plot(grids.new, sqrt(va))

#Estimated mean curve:
muest <- result_pir_acb$fitted_mean

# If we wanted to exclude players with only one measurement:
#more_than_one <- df_pir_acb %>% 
#  count(Name) %>%
#  filter(n > 1) %>%
#  select(Name) 

df_pir_acb1 <- df_pir_acb #%>% 
  #filter(Name %in% more_than_one$Name) 

L1 <- min(as.numeric(df_pir_acb1$Age)) # Range of the time.             
L2 <- max(as.numeric(df_pir_acb1$Age))

df_pir_acb2 <- df_pir_acb1 %>% 
  mutate(temp = (Age - L1) / (L2 - L1)) %>%
  mutate(temp = ifelse(temp < 0.00001, 0.00001, temp)) %>%
  mutate(temp = ifelse(temp > 0.99999, 0.99999, temp)) %>%
  select(Name, PIR, temp) %>%
  rename(Age = temp)

player_names_pir_acb <- as.character(unique(df_pir_acb2$Name))
#save(player_names_pir_acb, file = "player_names_pir_acb.RData")
#load("player_names_pir_acb.RData")

fpcs <- fpca.score(as.matrix(df_pir_acb2), grids.new, muest, evalest, eigenfest, sig2est, r)
pred <- fpca.pred(fpcs, muest, eigenfest) 
# Save table (fpcs will be used with FAA and FADA):
#write.table(fpcs, file = "fpcs.txt", col.names = FALSE, row.names = FALSE)
#fpcs <- read.table("Predictions/Files_from_diff_places/fpcs.txt")

# ----

# 3. FAA and FADA
library(Anthropometry)
X <- fpcs
numArch <- 10 
numRep <- 20
set.seed(2018)
s0 <- Sys.time()
lass_pir_acb <- stepArchetypesRawData(data = X, numArch = 1:numArch, numRep = numRep, verbose = FALSE)
s1 <- Sys.time() - s0
s1
#Time difference of 14.980789 mins
screeplot(lass_pir_acb)
#save(lass_pir_acb, file = "lass_pir_acb.RData")
#load("lass.RData")

numArchoid <- 4 # Number of archetypoids according to the screeplot.
s0 <- Sys.time()
res_ada_pir_acb <- archetypoids(numArchoid, data = X, step = FALSE, 
                           ArchObj = lass_pir_acb, nearest = "cand_ns", sequ = TRUE)
s1 <- Sys.time() - s0
s1
#Time difference of 32.729508 mins
#save(res_ada_pir_acb, file = "res_ada_pir_acb.RData")
#load("res_ada_pir_acb.RData")

str(res_ada_pir_acb)
# Archetypoids:
cases <- player_names_pir_acb[res_ada_pir_acb$cases]
cases
#[1] "Vazquez, Paco~AR4"   "Sabonis, A.~AAM"   "Bell, T.~B6W"   "Fernandez, Rudy~B0X" 
df_pir_acb1 %>%
  filter(Name %in% cases) %>%
  group_by(Name) %>%
  summarise(mean_pir = mean(PIR)) %>%
  arrange(desc(mean_pir))
#Sabonis, A.~AAM        28.4 
#Fernandez, Rudy~B0X    12.8 
#Vazquez, Paco~AR4      3.54
#Bell, T.~B6W           -7.50

#mean(df_pir_acb1[df_pir_acb1$Name == cases[2], "PIR"])

# Alphas:
dimnames(res_ada_pir_acb$alphas) <- list(player_names_pir_acb[res_ada_pir_acb$cases], player_names_pir_acb) 
alphas <- res_ada_pir_acb$alphas
names(which(alphas[2,] > 0.85)) 
#[1] "Berry, Walter~AAG"  "Magee, Kevin~AE7"  "Sabonis, A.~AAM" 

alphas[, grep("Gasol, Pau", colnames(alphas))]
#  Vazquez, Paco~AR4     Sabonis, A.~AAM        Bell, T.~B6W Fernandez, Rudy~B0X 
#         0.12753260          0.24944194          0.09914133          0.52387670
alphas[, grep("Rubio, Ricky", colnames(alphas))]
#Vazquez, Paco~AR4     Sabonis, A.~AAM        Bell, T.~B6W Fernandez, Rudy~B0X 
#        0.0000000           0.5985313           0.0000000           0.4099059
alphas[, grep("Biyombo", colnames(alphas))]
#Vazquez, Paco~AR4     Sabonis, A.~AAM        Bell, T.~B6W Fernandez, Rudy~B0X 
#       0.01305877          0.31239826          0.31169717          0.36278538
alphas[, grep("Porzingis", colnames(alphas))]
#  Vazquez, Paco~AR4     Sabonis, A.~AAM        Bell, T.~B6W Fernandez, Rudy~B0X 
#         0.03259843          0.25471587          0.26271462          0.44992671
alphas[, grep("Hezonja", colnames(alphas))]
# Vazquez, Paco~AR4     Sabonis, A.~AAM        Bell, T.~B6W Fernandez, Rudy~B0X 
#         0.1237563           0.2086463           0.3955814           0.2719861 
alphas[, grep("Doncic, Luka", colnames(alphas))]
# Vazquez, Paco~AR4     Sabonis, A.~AAM        Bell, T.~B6W Fernandez, Rudy~B0X 
#         0.0000000           0.5437567           0.0000000           0.4655219


# SET OF SIMILAR PLAYERS:
#colnames(alphas)[grep("Doncic, Luka", colnames(alphas))]
#[1] "Doncic, Luka~460"
player <- "Doncic, Luka~460"
source(".../Code_and_data_v1/do_similar_players.R")
do_similar_players(player, alphas, numArchoid, FALSE)
