library(dplyr)
library(tidyr)
library(ggplot2)

distance <-  as.matrix(read.csv("../../data/input/northeastamerica/distanceGridNEA.csv", row.names = 1, stringsAsFactors = FALSE))
#distance <-  as.matrix(read.csv("../../data/input/atlanticarea/distanceGridAA.csv", row.names = 1, stringsAsFactors = FALSE))

distance <- distance %>%
  replace(., col(.) == row(.), NA) %>% 
  as.data.frame() %>% mutate(destination = row.names(.)) %>% 
  pivot_longer(cols = -destination, names_to = 'departure', values_to = 'distance')


# true values 
meanInterDistance <- mean(distance$distance, na.rm = TRUE )
standardDeviationInterDistance <- sd(distance$distance, na.rm = TRUE )

# as in XML from Rougier's application
meanInterDistance = 300
standardDeviationInterDistance = 978

alpha0 = -2.9
alpha1 =  19.7


distance <- distance %>% mutate(logitW =  alpha0 - alpha1 * (distance - meanInterDistance) / standardDeviationInterDistance) %>% 
  mutate(W = 1/(1 + exp(-logitW))) 

  #filter(departure =='Pearl', destination == 'Escambia') 
  
distance %>% group_by(departure) %>% summarise(sumW = sum(W, na.rm = TRUE)) %>%
  #filter(departure == 'Potomac')
  ggplot(aes(x = sumW)) + geom_histogram() + 
  geom_vline(aes(xintercept = mean(sumW)), color="blue", linetype="dashed", size=1)
  

distance %>% group_by(departure) %>% summarise(sumW = sum(W, na.rm = TRUE)) %>%
  summarise(sumW = mean(sumW))