## Entering flow centrality simulations
## Files are in "/Users/Monnie/Dropbox/2018 Fall Leave/FlowThroughCentrality/Simulations"
## Data are in https://smu.app.box.com/folder/47024438738
setwd("../Nodes150Deg6Pert3")
library(tidyverse)
library(xtable)

# Enter all data files
varNames <- c("Network",paste0(c("Cent","Tier"), rep(1:150,each=2)))
dataBW <- read.csv(file="Nodes150_AvgDeg06_UW_01_Pert03_Betweennesss.csv",skip=5,nrows=100,header=FALSE)
colnames(dataBW) <- varNames
dataBW$Measure <- rep("BW",100)
dataCL <- read.csv(file="Nodes150_AvgDeg06_UW_01_Pert03_Closeness.csv",skip=5,nrows=100,header=FALSE)
colnames(dataCL) <- varNames
dataCL$Measure <- rep("CL",100)
dataFB <- read.csv(file="Nodes150_AvgDeg06_UW_01_Pert03_FlowBetween.csv",skip=5,nrows=100,header=FALSE)
colnames(dataFB) <- varNames
dataFB$Measure <- rep("FB",100)
dataFT <- read.csv(file="Nodes150_AvgDeg06_UW_01_Pert03_Flowthrough.csv",skip=5,nrows=100,header=FALSE)
colnames(dataFT) <- varNames
dataFT$Measure <- rep("FT",100)
dataSB <- read.csv(file="Nodes150_AvgDeg06_UW_01_Pert03_StableBetween.csv",skip=5,nrows=100,header=FALSE)
colnames(dataSB) <- varNames
dataSB$Measure <- rep("SB",100)

## Create tibbles with first six columns of each file
dataBWfirst6 <- dplyr::select(dataBW,c(Network:Tier6,Measure))
dataCLfirst6 <- dplyr::select(dataCL,c(Network:Tier6,Measure))
dataFBfirst6 <- dplyr::select(dataFB,c(Network:Tier6,Measure))
dataFTfirst6 <- dplyr::select(dataFT,c(Network:Tier6,Measure))
dataSBfirst6 <- dplyr::select(dataSB,c(Network:Tier6,Measure))
dataAllFirst6 <- bind_rows(dataBWfirst6,dataCLfirst6,dataFBfirst6,dataFTfirst6, dataSBfirst6)
dataAllFirst6$Measure <- as.factor(dataAllFirst6$Measure)

## Data in long format
tier6n150deg6pert3a <- gather(dataAllFirst6,key=CentLevel,value=Centrality,Cent1,Cent2,Cent3, Cent4, Cent5, Cent6,factor_key=TRUE) 
tier6n150deg6pert3b <- gather(dataAllFirst6,key=Level,value=Tier, Tier1,Tier2,Tier3,Tier4, Tier5, Tier6, factor_key=TRUE)
tier6n150deg6pert3 <-merge(tier6n150deg6pert3a,tier6n150deg6pert3b,by=c("Network","Measure"))
tier6n150deg6pert3 <- dplyr::select(tier6n150deg6pert3,-c(Tier1:Tier6,Cent1:Cent6))
tier6n150deg6pert3 <- separate(tier6n150deg6pert3,Network,into=c("Trash","Sample"),sep=" ")
tier6n150deg6pert3$Sample <- as.integer(tier6n150deg6pert3$Sample)
tier6n150deg6pert3 <- tier6n150deg6pert3[,-1]

pdf(file="tier6n150deg6pert3Lines.pdf")
ggplot(tier6n150deg6pert3,aes(x=Sample, y=Centrality, colour=CentLevel)) + geom_line() + facet_grid(Measure~.,scales="free")
dev.off()

xtable(tier6n150deg6pert3 %>% group_by(CentLevel,Measure) %>% summarise(mean=mean(Centrality),sd=sd(Centrality),cv=100*(sd(Centrality)/mean(Centrality))))

## More plots
http://www.sthda.com/english/wiki/ggplot2-error-bars-quick-start-guide-r-software-and-data-visualization

## Some histograms. Data are not Normal for the most part.
par(mfrow=c(5,3))
hist(dataAllFirst3$BWTier1)
hist(dataAllFirst3$BWTier2)
hist(dataAllFirst3$BWTier3)
hist(dataAllFirst3$FBTier1)
hist(dataAllFirst3$FBTier2)
hist(dataAllFirst3$FBTier3)
hist(dataAllFirst3$FTTier1)
hist(dataAllFirst3$FTTier2)
hist(dataAllFirst3$FTTier3)
hist(dataAllFirst3$SBTier1)
hist(dataAllFirst3$SBTier2)
hist(dataAllFirst3$SBTier3)
hist(dataAllFirst3$CLTier1)
hist(dataAllFirst3$CLTier2)
hist(dataAllFirst3$CLTier3)

## Create tibble in long format

## Run tests for variance on dependent variables

## References
## http://psychologicalstatistics.blogspot.com/2006/05/what-is-all-this-stuff-about.html
### https://www.rdocumentation.org/packages/vegan/versions/2.4-2/topics/mrpp
## http://cc.oulu.fi/~jarioksa/softhelp/vegan/html/mrpp.html
