#RCODE FOR GENERATING GRAPHS OF Pointing to the Ahmarian. Lithic technology and the el-Wad points of Al-Ansab 1 GENNAI, J., SCHEMMEL, M. & RICHTER, J.

#PRELIMINARY
install.packages("ggplot2")
library("ggplot2")
library("gridExtra")
library("ggsci")
install.packages("openxlsx")
library(openxlsx)
install.packages("dplyr")
install.packages('svglite')
library(dplyr)
library(scales)
install.packages("ggpubr")
library(ggpubr)
library(cowplot)
install.packages("RColorBrewer")
library(RColorBrewer)
library(forcats)
install.packages("rcartocolor")
library(rcartocolor)
install.packages("extrafont")
library(extrafont)
install.packages("ggstatsplot")
library(ggstatsplot)
install.packages("corrplot")
library(corrplot)
setwd("/Users/Jacopo/OneDrive/Ansab Paper")#C:\Users\Jacopo\OneDrive\Ansab Paper
getwd()


#DIMENSIONS BLANKS AND NEGATIVES 
AnsabBlanks <- read.csv(file.choose()) #take the Dimensions Ansab file
str(AnsabBlanks)
AnsabBlanks1 <- subset(AnsabBlanks, StructureCat %in% c("Blade", "Bladelet", "Flake", "Blade and Bladelet", "LCN Core", "LCN Debitage", "El-Wad"), select=c(Entirety, StructureCat, EI, Length, Width, Thick, Curv, TechCat))
#REORDER THE CATEGORIES 
AnsabBlanks1$TechCat_F <- factor(AnsabBlanks1$TechCat, levels = c(
"Simple", "mgmt", "Spall", "Tablet", "Init", "LCN Debitage", "LCN Core", "El-Wad"))
AnsabBlanks1$StructureCat_F <- factor(AnsabBlanks1$StructureCat, levels = c(
"Blade", "Bladelet", "El-Wad", "Blade and Bladelet", "Flake", "LCN Debitage", "LCN Core"))
str(AnsabBlanks1)
#LENGTH
AnsabBlanksLength <- subset(AnsabBlanks1, Entirety %in% c("Complete"), select=c(StructureCat, EI, Length, Width, Thick, Curv, TechCat))
AnsabBlanksLength$TechCat_F <- factor(AnsabBlanksLength$TechCat, levels = c(
"Simple", "mgmt", "Spall", "Tablet", "Init", "LCN Debitage", "LCN Core", "El-Wad"))
str(AnsabBlanksLength)
AnsabBlanksLength$StructureCat_F <- factor(AnsabBlanksLength$StructureCat, levels = c(
"Blade", "Bladelet", "El-Wad", "Blade and Bladelet", "Flake", "LCN Debitage", "LCN Core"))
str(AnsabBlanksLength)
#VIOLIN PLOT
Pp<-ggplot(AnsabBlanksLength, aes(TechCat_F, y=Length, fill=TechCat_F))+
geom_violin()+
theme_classic()+
xlab("")+
ylab("mm")+
facet_grid(~StructureCat_F, labeller = labeller(StructureCat_F=label_wrap_gen(width = 8)), scales="free", space="free")+
theme(legend.position="bottom") +
theme(
axis.text=element_text(face="bold", size=10, color="black", family="Arial"),
axis.text.y=element_text(face="bold", size=10, family="Arial"),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
legend.title=element_blank(),
legend.text=element_text(size=10, family="Arial"),
strip.text = element_text(face="bold", size=11, family="Arial"))+
scale_y_continuous(breaks= seq(0, 120, by = 10), limits = c(0, 120))+
scale_fill_manual(values=c("#999999", "#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00", "#CC79A7")) +
guides(fill=guide_legend(nrow=1,byrow=TRUE))
Pp
#ADDING THE BOXPLOT
Ppcombined <-Pp+geom_boxplot(width=0.2, colour="black")
Ppcombined
dev.size(units=c("cm"))
#WIDTH
#VIOLIN PLOT
Pp1<-ggplot(AnsabBlanks1, aes(TechCat_F, y=Width, fill=TechCat_F))+
geom_violin()+
theme_classic()+
xlab("")+
ylab("mm")+
facet_grid(~StructureCat_F, labeller = labeller(StructureCat_F=label_wrap_gen(width = 8)), scales = "free")+
theme(legend.position="bottom") +
theme(
axis.text=element_text(face="bold", size=10, color="black", family="Arial"),
axis.text.y=element_text(face="bold", size=10, family="Arial"),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
legend.title=element_blank(),
legend.text=element_text(size=10, family="Arial"),
strip.text = element_text(face="bold", size=10, family="Arial"))+
scale_y_continuous(breaks= seq(0, 90, by = 10), limits = c(0, 90))+
scale_fill_manual(values=c("#999999", "#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00", "#CC79A7")) +
geom_hline(yintercept=12, linetype="dashed", color = "black")+
guides(fill=guide_legend(nrow=1,byrow=TRUE))
Pp1
#ADDING THE BOXPLOT
Pp1combined <-Pp1+geom_boxplot(width=0.2, colour="black")
Pp1combined
dev.size()
#THICKNESS
#SELECT ONLY DEBITAGE CATEGORIES WITH TICKNESS VALUES
AnsabBlanks2 <- subset(AnsabBlanks1, StructureCat_F %in% c("Blade", "Bladelet", "Flake", "Blade and Bladelet", "El-Wad"), select=c(StructureCat_F, Thick, TechCat_F))
#VIOLIN PLOT
Pp2<-ggplot(AnsabBlanks2, aes(TechCat_F, y=Thick, fill=TechCat_F))+
geom_violin()+
theme_classic()+
xlab("")+
ylab("mm")+
facet_grid(~StructureCat_F, labeller = labeller(StructureCat_F=label_wrap_gen(width = 15)), scales = "free", space="free")+
theme(legend.position="bottom") +
theme(
axis.text=element_text(face="bold", size=10, color="black", family="Arial"),
axis.text.y=element_text(face="bold", size=10, family="Arial"),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
legend.title=element_blank(),
legend.text=element_text(size=10, family="Arial"),
strip.text = element_text(face="bold", size=10, family="Arial"))+
scale_y_continuous(breaks= seq(0, 50, by = 5), limits = c(0, 50))+
scale_fill_manual(values=c("#999999", "#E69F00", "#F0E442", "#0072B2", "#D55E00", "#CC79A7"))+
guides(fill=guide_legend(nrow=1,byrow=TRUE))
Pp2
#ADDING THE BOXPLOT
Pp2combined <-Pp2+geom_boxplot(width=0.2, colour="black")
Pp2combined
dev.size(units = c("cm"))
#EXPORTING
ggsave("AnsabLength.tiff", plot=Ppcombined, units="cm", width=50.80000, height=26.51125, dpi=600, compression = 'lzw')
ggsave("AnsabLength.svg", plot=Ppcombined, width=50.80000, height=26.51125, units = "cm")
ggsave("AnsabWidth.tiff", plot=Pp1combined, units="cm", width=50.80000, height=26.51125, dpi=600, compression = 'lzw')
ggsave("AnsabWidth.svg", plot=Pp1combined, width=50.80000, height=26.51125, units = "cm")
ggsave("AnsabThick.tiff", plot=Pp2combined, units="cm", width=50.80000, height=26.51125, dpi=600, compression = 'lzw')
ggsave("AnsabThick.svg", plot=Pp2combined, width=50.80000, height=26.51125, units = "cm")

#CORES 
AnsabCores <- subset(AnsabBlanks, StructureCat %in% c("Core"), select=c(Entirety, StructureCat, Length, Width, Thick, TechCat))
#SELECTING ONLY COMPLETE CORES
AnsabCores1 <- subset(AnsabCores, Entirety %in% c("VERO"), select=c(StructureCat, Length, Width, Thick, TechCat))
#REORDERING THE CORES CATEGORIES
AnsabCores1$TechCat_F <- factor(AnsabCores1$TechCat, levels = c(
"Semi Tournant", "Parallel Edges", "Trans. Carinated", "Narrow Fronted", "sur Tranche", "Pre-Core", "Others"))
str(AnsabCores1)
#VIOLIN PLOT
Ppcores<-ggplot(AnsabCores1, aes(TechCat_F, y=Length, fill=TechCat_F))+
geom_violin()+
theme_classic()+
xlab("")+
ylab("mm")+
facet_grid(~TechCat_F, labeller = labeller(TechCat_F=label_wrap_gen(width = 8)), scales="free", space="free")+
theme(legend.position="bottom") +
theme(
axis.text=element_text(face="bold", size=10, color="black", family="Arial"),
axis.text.y=element_text(face="bold", size=10, family="Arial"),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
legend.title=element_blank(),
legend.text=element_text(size=11, family="Times New Roman"),
strip.text = element_text(face="bold", size=10, family="Arial"))+
scale_y_continuous(breaks= seq(0, 120, by = 10), limits = c(0, 120))+
scale_fill_brewer(palette="Set2")+
guides(fill=guide_legend(nrow=1,byrow=TRUE))
Ppcores
#ADDING THE BOXPLOT
Ppcorescombined <-Ppcores+geom_boxplot(width=0.2, colour="black")
Ppcorescombined
#WIDTH
#VIOLIN PLOT
Pp1cores<-ggplot(AnsabCores1, aes(TechCat_F, y=Width, fill=TechCat_F))+
geom_violin()+
theme_classic()+
xlab("")+
ylab("mm")+
facet_grid(~TechCat_F, labeller = labeller(TechCat_F=label_wrap_gen(width = 8)), scales="free", space="free")+
theme(legend.position="bottom") +
theme(
axis.text=element_text(face="bold", size=10, color="black", family="Arial"),
axis.text.y=element_text(face="bold", size=10, family="Arial"),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
legend.title=element_blank(),
legend.text=element_text(size=11, family="Arial"),
strip.text = element_text(face="bold", size=10, family="Arial"))+
scale_y_continuous(breaks= seq(0, 100, by = 10), limits = c(0, 100))+
scale_fill_brewer(palette="Set2")+
guides(fill=guide_legend(nrow=1,byrow=TRUE))
Pp1cores
#ADDING THE BOXPLOT
Pp1corescombined <-Pp1cores+geom_boxplot(width=0.2, colour="black")
Pp1corescombined
#THICKNESS
#VIOLIN PLOT
Pp2cores<-ggplot(AnsabCores1, aes(TechCat_F, y=Thick, fill=TechCat_F))+
geom_violin()+
theme_classic()+
xlab("")+
ylab("mm")+
facet_grid(~TechCat_F, labeller = labeller(TechCat_F=label_wrap_gen(width = 8)), scales="free", space="free")+
theme(legend.position="bottom") +
theme(
axis.text=element_text(face="bold", size=10, color="black", family="Arial"),
axis.text.y=element_text(face="bold", size=10, family="Arial"),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
legend.title=element_blank(),
legend.text=element_text(size=10, family="Arial"),
strip.text = element_text(face="bold", size=10, family="Arial"))+
scale_y_continuous(breaks= seq(0, 120, by = 10), limits = c(0, 120))+
scale_fill_brewer(palette="Set2")+
guides(fill=guide_legend(nrow=1,byrow=TRUE))
Pp2cores
#ADDING THE BOXPLOT
Pp2corescombined <-Pp2cores+geom_boxplot(width=0.2, colour="black")
Pp2corescombined
dev.size(units=c("cm"))
#EXPORTING
ggsave("AnsabCoresLength.tiff", plot=Ppcorescombined, units="cm", width= 50.80000, height=26.51125, dpi=600, compression = 'lzw')
ggsave("AnsabCoresLength.svg", plot=Ppcorescombined, width= 50.80000, height=26.51125, units = "cm")
ggsave("AnsabCoresWidth.tiff", plot=Pp1corescombined, units="cm", width= 50.80000, height=26.51125,, dpi=600, compression = 'lzw')
ggsave("AnsabCoresWidth.svg", plot=Pp1corescombined, width= 50.80000, height=26.51125,, units = "cm")
ggsave("AnsabCoresThick.tiff", plot=Pp2corescombined, units="cm", width= 50.80000, height=26.51125,, dpi=600, compression = 'lzw')
ggsave("AnsabCoresThick.svg", plot=Pp2corescombined, width= 50.80000, height=26.51125, units = "cm")
#SUMMARY TABLES
dfDebitageL <- AnsabBlanksLength %>%
group_by(StructureCat_F, TechCat_F) %>%
  summarise(
    count = n(),
    MIN= min(Length, na.rm = TRUE),
    IQ = quantile(Length, 0.25, na.rm = TRUE),
    median = median(Length, na.rm = TRUE),
    IIIQ= quantile(Length, 0.75, na.rm = TRUE),
    MAX = max(Length, na.rm = TRUE))
print(tbl_df(dfDebitageL), n=24)
df1DebitageL <- as.data.frame(dfDebitageL)
df1DebitageL
write.xlsx(dfDebitageL, "AnsabDebitage_Length.xlsx")
df1DebitageW <- AnsabBlanks1 %>%
group_by(StructureCat_F, TechCat_F) %>%
  summarise(
    count = n(),
    MIN= min(Width, na.rm = TRUE),
    IQ = quantile(Width, 0.25, na.rm = TRUE),
    median = median(Width, na.rm = TRUE),
    IIIQ= quantile(Width, 0.75, na.rm = TRUE),
    MAX = max(Width, na.rm = TRUE))
print(tbl_df(df1DebitageW), n=24)
df1DebitagesW <- as.data.frame(dfDebitageW)
df1DebitageW
write.xlsx(df1DebitageW, "AnsabDebitage_Width.xlsx")
df1DebitageT <- AnsabBlanks2 %>%
group_by(StructureCat_F, TechCat_F) %>%
  summarise(
    count = n(),
    MIN= min(Thick, na.rm = TRUE),
    IQ = quantile(Thick, 0.25, na.rm = TRUE),
    median = median(Thick, na.rm = TRUE),
    IIIQ= quantile(Thick, 0.75, na.rm = TRUE),
    MAX = max(Thick, na.rm = TRUE))
print(tbl_df(dfDebitageT), n=24)
df1DebitageT <- as.data.frame(df1DebitageT)
df1DebitageT
write.xlsx(df1DebitageT, "AnsabDebitage_Thick.xlsx")

dfCoresL <- 
group_by(AnsabCores1, TechCat_F) %>%
  summarise(
    count = n(),
    MIN= min(Length, na.rm = TRUE),
    IQ = quantile(Length, 0.25, na.rm = TRUE),
    median = median(Length, na.rm = TRUE),
    IIIQ= quantile(Length, 0.75, na.rm = TRUE),
    MAX = max(Length, na.rm = TRUE))
print(tbl_df(dfCoresL), n=24)
dfCoresL <- as.data.frame(dfCoresL)
dfCoresL
write.xlsx(dfCoresL, "AnsabCores_Length.xlsx")
dfCoresW <- 
group_by(AnsabCores1, TechCat_F) %>%
  summarise(
    count = n(),
    MIN= min(Width, na.rm = TRUE),
    IQ = quantile(Width, 0.25, na.rm = TRUE),
    median = median(Width, na.rm = TRUE),
    IIIQ= quantile(Width, 0.75, na.rm = TRUE),
    MAX = max(Width, na.rm = TRUE))
print(tbl_df(dfCoresW), n=24)
dfCoresW <- as.data.frame(dfCoresW)
dfCoresW
write.xlsx(dfCoresW, "AnsabCores_Width.xlsx")
dfCoresT <- 
group_by(AnsabCores1, TechCat_F) %>%
  summarise(
    count = n(),
    MIN= min(Thick, na.rm = TRUE),
    IQ = quantile(Thick, 0.25, na.rm = TRUE),
    median = median(Thick, na.rm = TRUE),
    IIIQ= quantile(Thick, 0.75, na.rm = TRUE),
    MAX = max(Thick, na.rm = TRUE))
print(tbl_df(dfCoresT), n=24)
dfCoresT <- as.data.frame(dfCoresT)
dfCoresT
write.xlsx(dfCoresT, "AnsabCores_Thick.xlsx")

#EI
str(AnsabBlanksLength)
#VIOLIN PLOT
PpEI<-ggplot(AnsabBlanksLength, aes(TechCat_F, y=EI, fill=TechCat_F))+
geom_violin()+
theme_classic()+
xlab("")+
ylab("Elongation Index")+
facet_grid(~StructureCat_F, labeller = labeller(StructureCat_F=label_wrap_gen(width = 8)), scales="free", space="free")+
theme(legend.position="bottom") +
theme(
axis.text=element_text(face="bold", size=10, color="black", family="Arial"),
axis.text.y=element_text(face="bold", size=10, family="Arial"),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
legend.title=element_blank(),
legend.text=element_text(size=10, family="Arial"),
strip.text = element_text(face="bold", size=10, family="Arial"))+
scale_y_continuous(breaks= seq(0, 12, by = 1), limits = c(0, 12))+
scale_fill_manual(values=c("#999999", "#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00", "#CC79A7"))+
guides(fill=guide_legend(nrow=1,byrow=TRUE))
PpEI
#ADDING BOXPLOT
PpEIcombined <-PpEI+geom_boxplot(width=0.1, colour="black")
PpEIcombined
dev.size(units=c("cm"))
#EXPORTING
ggsave("AnsabElongation.tiff", plot=PpEIcombined, units="cm", width=50.80000, height=26.51125, dpi=600, compression = 'lzw')
ggsave("AnsabElongation.svg", plot=PpEIcombined, units="cm", width=50.80000, height=26.51125)

df1EI <- AnsabBlanksLength %>%
group_by(StructureCat_F, TechCat_F) %>%
  summarise(
    count = n(),
    MIN= min(EI, na.rm = TRUE),
    IQ = quantile(EI, 0.25, na.rm = TRUE),
    median = median(EI, na.rm = TRUE),
    IIIQ= quantile(EI, 0.75, na.rm = TRUE),
    MAX = max(EI, na.rm = TRUE))
print(tbl_df(df1EI), n=24)
df1EI <- as.data.frame(df1EI)
df1EI
write.xlsx(df1EI, "AnsabDebitage_EI.xlsx")


#CURVATURE
AnsabCurvature <- subset(AnsabBlanksLength, StructureCat %in% c("Blade", "Bladelet", "Blade and Bladelet", "El-Wad"), select=c(StructureCat, EI, Curv, Length, Width, Thick, TechCat))
str(AnsabCurvature)
AnsabCurvature
#REORDERING CATEGORIES
AnsabCurvature$TechCat_F <- factor(AnsabCurvature$TechCat, levels = c(
"Simple", "mgmt", "Spall", "El-Wad"))
AnsabCurvature$StructureCat_F <- factor(AnsabCurvature$StructureCat, levels = c(
"Blade", "Bladelet", "Blade and Bladelet", "El-Wad"))
str(AnsabCurvature)
pal <- (values=c("#999999", "#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00", "#CC79A7"))
show_col(pal)

#BOXPLOT
Curv<-ggplot(AnsabCurvature, aes(StructureCat_F, y=Curv, fill=TechCat_F))+
geom_boxplot()+
theme_classic()+
xlab("")+
ylab("Curvature")+
theme(legend.position="bottom") +
facet_grid(~StructureCat_F, labeller = labeller(StructureCat_F=label_wrap_gen(width = 8)), scales = "free", space="free")+
theme(
axis.text=element_text(face="bold", size=10, color="black", family="Arial"),
axis.text.y=element_text(face="bold", size=10, family="Arial"),
axis.text.x=element_blank(),
axis.ticks.x=element_blank(),
legend.title=element_blank(),
legend.text=element_text(size=11, family="Arial"),
strip.text = element_text(face="bold", size=11, family="Arial"))+
scale_y_continuous(breaks= c(0, 2.9, 5.9, 8.9, 14.9, 20), limits = c(0, 22))+
scale_fill_manual(values=c("#999999", "#E69F00", "#56B4E9", "#CC79A7"))+ 
guides(fill=guide_legend(nrow=1,byrow=TRUE))+
geom_hline(yintercept=2.9, linetype="dashed", color = "black")+
geom_hline(yintercept=5.9, linetype="dashed", color = "black")+
geom_hline(yintercept=8.9, linetype="dashed", color = "black")+
geom_hline(yintercept=14.9, linetype="dashed", color = "black")
Curv
dev.size(units="cm")
#EXPORTING
ggsave("AnsabCurvature.tiff", plot=Curv, units="cm", width=17.78000, height=17.75354, dpi=600, compression = 'lzw')
ggsave("AnsabCurvature.svg", plot=Curv, units="cm", width=17.78000, height=17.75354)

df1Curv <- AnsabCurvature %>%
group_by(StructureCat_F, TechCat_F) %>%
  summarise(
    count = n(),
    MIN= min(Curv, na.rm = TRUE),
    IQ = quantile(Curv, 0.25, na.rm=TRUE),
    median = median(Curv, na.rm = TRUE),
    IIIQ= quantile(Curv, 0.75, na.rm=TRUE),
    MAX = max(Curv, na.rm = TRUE))
print(tbl_df(dfCurv), n=24)
df1Curv <- as.data.frame(df1Curv)
df1Curv
write.xlsx(df1Curv, "AnsabDebitage_Curv.xlsx")

#WILCOXON TEST UNPAIRED BETWEEN SIMPLE BLADES AND BLADELETS AND EL-WAD POINTS
#Selecting the simple complete blanks
Simple_ElWad_Complete<- subset(AnsabCurvature, TechCat %in% c("Simple", "El-Wad"), select=c(StructureCat, Curv, EI, Length, Width, Thick, TechCat))
str(Simple_ElWad_Complete)
#selecting simple complete and semi blanks
BladesBladeletsBlanks <- subset(AnsabBlanks, StructureCat %in% c("Blade", "Bladelet", "Blade and Bladelet", "El-Wad"), select=c(StructureCat, Curv, Length, Width, Thick, TechCat))
Simple_ElWad_CompleteSemi<- subset(BladesBladeletsBlanks, TechCat %in% c("Simple", "El-Wad"), select=c(StructureCat, Curv, Length, Width, Thick, TechCat))
str(Simple_ElWad_CompleteSemi)
install.packages("gginnards")
library(gginnards)

LengthBladeletsBladesElWad <- ggstatsplot::ggbetweenstats(
  data = filter(Simple_ElWad_CompleteSemi, StructureCat %in% c("Blade and Bladelet", "El-Wad")),
  x = StructureCat,
  y = Length,
  violin.args = list(width = 0.5, alpha = 0),
  centrality.plotting = TRUE, # whether mean for each group is to be displayed
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "mm", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)
LengthBladeletsBladesElWad 

## delete the geom_point layer
LengthBladeletsBladesElWad <- gginnards::delete_layers(x = LengthBladeletsBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
LengthBladeletsBladesElWad1 <- LengthBladeletsBladesElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#9EBCDA", "#000000"))
LengthBladeletsBladesElWad1


WidthBladeletsBladesElWad <- ggstatsplot::ggbetweenstats(
  data = filter(Simple_ElWad_CompleteSemi, StructureCat %in% c("Blade and Bladelet", "El-Wad")),
  x = StructureCat,
  y = Width,
  violin.args = list(width = 0.5, alpha = 0),
  centrality.plotting = TRUE, # whether mean for each group is to be displayed
  medianan.label.size = 2, 
  type = "nonparametric",  
  k = 1, 
  xlab = "",
  ylab = "mm", 
  ggtheme = ggplot2::theme_classic(),
 ggplot.component = ggplot2::theme(axis.ticks.x=element_blank(), 
					     axis.text=element_text(face="bold", size=10, color="black", family="Arial"))
)
WidthBladeletsBladesElWad

## delete the geom_point layer
WidthBladeletsBladesElWad<- gginnards::delete_layers(x = WidthBladeletsBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
WidthBladeletsBladesElWad1 <- WidthBladeletsBladesElWad+ 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#9EBCDA", "#000000"))
WidthBladeletsBladesElWad1
 


ThickBladeletsBladesElWad <- ggstatsplot::ggbetweenstats(
  data = filter(Simple_ElWad_CompleteSemi, StructureCat %in% c("Blade and Bladelet", "El-Wad")),
  x = StructureCat,
  y = Thick,
point.args = list(position = ggplot2::position_jitterdodge(dodge.width = 0.75), alpha
    = 1, size = 1, stroke = 1),
  violin.args = list(width = 0.5, alpha = 0),
  centrality.plotting = TRUE, # whether mean for each group is to be displayed
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "mm", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)+
ggplot2::scale_color_manual(values =c("#9EBCDA", "#000000"))
ThickBladeletsBladesElWad

## delete the geom_point layer
ThickBladeletsBladesElWad<- gginnards::delete_layers(x = ThickBladeletsBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
ThickBladeletsBladesElWad1 <- ThickBladeletsBladesElWad+ 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#9EBCDA", "#000000"))
ThickBladeletsBladesElWad1


EIBladeletsBladesElWad <- ggstatsplot::ggbetweenstats(
  data = filter(Simple_ElWad_Complete, StructureCat %in% c("Blade and Bladelet", "El-Wad")),
  x = StructureCat,
  y = EI,
point.args = list(position = ggplot2::position_jitterdodge(dodge.width = 0.75), alpha
    = 1, size = 1, stroke = 1),
  violin.args = list(width = 0.5, alpha = 0),
  centrality.plotting = TRUE, # whether mean for each group is to be displayed
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "Elongation Index", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)+
ggplot2::scale_color_manual(values =c("#9EBCDA", "#000000"))
EIBladeletsBladesElWad

EIBladeletsBladesElWad<- gginnards::delete_layers(x = EIBladeletsBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
EIBladeletsBladesElWad1 <- EIBladeletsBladesElWad+ 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#9EBCDA", "#000000"))
EIBladeletsBladesElWad1



CurvBladeletsBladesElWad <- ggstatsplot::ggbetweenstats(
  data = filter(Simple_ElWad_Complete, StructureCat %in% c("Blade and Bladelet", "El-Wad")),
  x = StructureCat,
  y = Curv,
point.args = list(position = ggplot2::position_jitterdodge(dodge.width = 0.75), alpha
    = 1, size = 1, stroke = 1),
  violin.args = list(width = 0.5, alpha = 0),
  centrality.plotting = TRUE, # whether mean for each group is to be displayed
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "Curvature", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)+
ggplot2::scale_color_manual(values =c("#9EBCDA", "#000000"))
CurvBladeletsBladesElWad

CurvBladeletsBladesElWad<- gginnards::delete_layers(x = CurvBladeletsBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
CurvBladeletsBladesElWad1 <- CurvBladeletsBladesElWad+ 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#9EBCDA", "#000000"))
CurvBladeletsBladesElWad1


BladeBladeletElwadarranged <- ggarrange(LengthBladeletsBladesElWad1, 
WidthBladeletsBladesElWad1, 
ThickBladeletsBladesElWad1,
CurvBladeletsBladesElWad1,
EIBladeletsBladesElWad1, ncol = 3, nrow = 2, labels  = "AUTO")
BladeBladeletElwadarranged
dev.size(unit=c("cm"))
ggsave("Wilcoxon_BladeBladeletElwad.svg", plot=BladeBladeletElwadarranged, units="cm", width=50.80000, height=26.51125)
ggsave("Wilcoxon_BladeBladeletElwad.tiff", plot=BladeBladeletElwadarranged, units="cm", width=50.80000, height=26.51125, dpi=600, compression = 'lzw')


#SIMPLE BLADES AND EL-WAD WILCOXON TEST
LengthBladesElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_Complete, StructureCat %in% c("Blade", "El-Wad")),
  x = StructureCat,
  y = Length,
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "mm", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = list(ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial"))))
LengthBladesElWad 

## delete the geom_point layer
LengthBladesElWad <- gginnards::delete_layers(x = LengthBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
LengthBladesElWad1 <- LengthBladesElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#FDB863", "#000000"))
LengthBladesElWad1

WidthBladesElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_CompleteSemi, StructureCat %in% c("Blade", "El-Wad")),
  x = StructureCat,
  y = Width,
  violin.args = list(width = 0.5, alpha = 0),
  median.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "mm", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
) 
WidthBladesElWad 

## delete the geom_point layer
WidthBladesElWad <- gginnards::delete_layers(x = WidthBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
WidthBladesElWad1 <- WidthBladesElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#FDB863", "#000000"))
WidthBladesElWad1

ThickBladesElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_CompleteSemi, StructureCat %in% c("Blade", "El-Wad")),
  x = StructureCat,
  y = Thick,
  violin.args = list(width = 0.5, alpha = 0),
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "mm", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)
ThickBladesElWad
## delete the geom_point layer
ThickBladesElWad <- gginnards::delete_layers(x = ThickBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
ThickBladesElWad1 <- ThickBladesElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#FDB863", "#000000"))
ThickBladesElWad1

EIBladesElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_Complete, StructureCat %in% c("Blade", "El-Wad")),
  x = StructureCat,
  y = EI,
  violin.args = list(width = 0.5, alpha = 0),
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "Elongation", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)
EIBladesElWad
## delete the geom_point layer
EIBladesElWad <- gginnards::delete_layers(x = EIBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
EIBladesElWad1 <- EIBladesElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#FDB863", "#000000"))
EIBladesElWad1

CurvBladesElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_Complete, StructureCat %in% c("Blade", "El-Wad")),
  x = StructureCat,
  y = Curv,
  violin.args = list(width = 0.5, alpha = 0),
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "Curvature", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)
CurvBladesElWad
## delete the geom_point layer
CurvBladesElWad <- gginnards::delete_layers(x = CurvBladesElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
CurvBladesElWad1 <- CurvBladesElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#FDB863", "#000000"))
CurvBladesElWad1


BladeElwadarranged <- ggarrange(LengthBladesElWad1, 
WidthBladesElWad1, 
ThickBladesElWad1,
CurvBladesElWad1,
EIBladesElWad1, ncol = 3, nrow = 2, labels  = "AUTO")
BladeElwadarranged
dev.size(unit=c("cm"))
ggsave("Wilcoxon_BladeElwad.svg", plot=BladeElwadarranged, units="cm", width=50.80000, height=26.51125)
ggsave("Wilcoxon_BladeElwad.tiff", plot=BladeElwadarranged, units="cm", width=50.80000, height=26.51125, dpi=600, compression = 'lzw')


#SIMPLE BLADELETS AND EL-WAD POINTS WILCOXON TEST
LengthBladeletsElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_Complete, StructureCat %in% c("Bladelet", "El-Wad")),
  x = StructureCat,
  y = Length,
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "mm", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = list(ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial"))))
LengthBladeletsElWad 

## delete the geom_point layer
LengthBladeletsElWad <- gginnards::delete_layers(x = LengthBladeletsElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
LengthBladeletsElWad1 <- LengthBladeletsElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#00aa00ff", "#000000"))
LengthBladeletsElWad1

WidthBladeletsElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_CompleteSemi, StructureCat %in% c("Bladelet", "El-Wad")),
  x = StructureCat,
  y = Width,
  violin.args = list(width = 0.5, alpha = 0),
  median.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "mm", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
) 
WidthBladeletsElWad 

## delete the geom_point layer
WidthBladeletsElWad <- gginnards::delete_layers(x = WidthBladeletsElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
WidthBladeletsElWad1 <- WidthBladeletsElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#00aa00ff", "#000000"))
WidthBladeletsElWad1

ThickBladeletsElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_CompleteSemi, StructureCat %in% c("Bladelet", "El-Wad")),
  x = StructureCat,
  y = Thick,
  violin.args = list(width = 0.5, alpha = 0),
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "mm", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)
ThickBladeletsElWad
## delete the geom_point layer
ThickBladeletsElWad <- gginnards::delete_layers(x = ThickBladeletsElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
ThickBladeletsElWad1 <- ThickBladeletsElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#00aa00ff", "#000000"))
ThickBladeletsElWad1

EIBladeletsElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_Complete, StructureCat %in% c("Bladelet", "El-Wad")),
  x = StructureCat,
  y = EI,
  violin.args = list(width = 0.5, alpha = 0),
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "Elongation", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)
EIBladeletsElWad
## delete the geom_point layer
EIBladeletsElWad <- gginnards::delete_layers(x = EIBladeletsElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
EIBladeletsElWad1 <- EIBladeletsElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#00aa00ff", "#000000"))
EIBladeletsElWad1

CurvBladeletsElWad <- ggstatsplot::ggbetweenstats(
  data= filter(Simple_ElWad_Complete, StructureCat %in% c("Bladelet", "El-Wad")),
  x = StructureCat,
  y = Curv,
  violin.args = list(width = 0.5, alpha = 0),
  medianan.label.size = 2, # size of the label for mean
  type = "nonparametric", # which type of test is to be run
  k = 1, # number of decimal places for statistical results
  xlab = "",
  ylab = "Curvature", # label for the y-axis variable
   ggtheme = ggplot2::theme_classic(),# choosing a different theme
 ggplot.component = ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
)
CurvBladeletsElWad
## delete the geom_point layer
CurvBladeletsElWad <- gginnards::delete_layers(x = CurvBladeletsElWad, match_type = "GeomPoint")
## add a new layers for points with a different shape and colour
CurvBladeletsElWad1 <- CurvBladeletsElWad + 
geom_point(shape = 21, size=2, position = position_jitterdodge(), aes(color = StructureCat, fill=StructureCat))+
scale_color_manual(values=c("#000000","#000000"))+
scale_fill_manual(values =c("#00aa00ff", "#000000"))
CurvBladeletsElWad1


BladeletElwadarranged <- ggarrange(LengthBladeletsElWad1, 
WidthBladeletsElWad1, 
ThickBladeletsElWad1,
CurvBladeletsElWad1,
EIBladeletsElWad1, ncol = 3, nrow = 2, labels  = "AUTO")
BladeletElwadarranged
dev.size(unit=c("cm"))
ggsave("Wilcoxon_BladeletElwad.svg", plot=BladeletElwadarranged, units="cm", width=50.80000, height=26.51125)
ggsave("Wilcoxon_BladeletElwad.tiff", plot=BladeletElwadarranged, units="cm", width=50.80000, height=26.51125, dpi=600, compression = 'lzw')


#MORPHOTECHNICAL ATTRIBUTES TEST
MorphoAttributes <- read.csv(file.choose()) #take the morphotechnological attributes Ansab
str(MorphoAttributes)
MorphoAttributes$Promod_F <- factor(MorphoAttributes$Promod, levels = c(
"Straight and Slightly Curved", "Curved and Very Curved", "Twisted"))
str(MorphoAttributes)

display.brewer.pal(11, "RdYlBu")
brewer.pal(11, "RdYlBu")

plot1 <-ggbarstats(
  data             = MorphoAttributes,
  x                = Promod_F,
  y                = StructureCat,
  type="p",
bf.message=FALSE,
xlab = "",
ylab = "Profile",
legend.title="Profile",
ggtheme = ggplot2::theme_classic(),# choosing a different theme
ggplot.component = list(
ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
ggplot2::scale_fill_manual(values =c("#A50026", "#FFFFBF", "#4575B4")))
)
plot1

plot2 <-ggbarstats(
data= MorphoAttributes,
  x= Shape,
  y= StructureCat,
  type="p",
bf.message=FALSE,
xlab = "",
ylab = "Shape",
legend.title="Shape",
ggtheme = ggplot2::theme_classic(),# choosing a different theme
ggplot.component = list(
ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
ggplot2::scale_fill_manual(values =c("#A50026", "#F46D43", "#FFFFBF", "#ABD9E9", "#4575B4")))
)
plot2

plot3 <-
ggbarstats(
  data= MorphoAttributes,
  x= Symmetry,
  y= StructureCat,
  type="p",
bf.message=FALSE,
xlab = "",
ylab = "Symmetry",
legend.title="Symmetry",
ggtheme = ggplot2::theme_classic(),# choosing a different theme
ggplot.component = list(
ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
ggplot2::scale_fill_manual(values =c("#A50026", "#F46D43", "#FFFFBF", "#ABD9E9", "#4575B4")))
)
plot3

plot4 <-
ggbarstats(
  data             = filter(MorphoAttributes, DEndMorph %in% c("Feathered", "Plunging", "Hinged", "Stepped", "ND")),
  x                = DEndMorph,
  y                = StructureCat,
  type="p",
bf.message=FALSE,
xlab = "",
  ylab = "Distal end",
legend.title="Distal end",
ggtheme = ggplot2::theme_classic(),# choosing a different theme
ggplot.component = list(
ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
ggplot2::scale_fill_manual(values =c("#A50026", "#F46D43", "#FFFFBF", "#ABD9E9", "#4575B4")))
)
plot4

plot5 <-
ggbarstats(
  data             = filter(MorphoAttributes, Out %in% c("(sub)Parallel Edges", "Convergent", "Off-axis", "ND")),
  x                = Out,
  y                = StructureCat,

  type="p",
bf.message=FALSE,
xlab = "",
  ylab = "Outline",
legend.title="Outline",
ggtheme = ggplot2::theme_classic(),# choosing a different theme
ggplot.component = list(
ggplot2::theme( 
axis.ticks.x=element_blank(), 
axis.text=element_text(face="bold", size=10, color="black", family="Arial")),
ggplot2::scale_fill_manual(values =c("#A50026", "#F46D43", "#ABD9E9", "#4575B4")))
)
plot5

Morphochisquarearranged <- ggarrange(plot1, 
plot2, 
plot3,
plot4,
plot5, ncol = 3, nrow = 2, labels  = "AUTO")
Morphochisquarearranged
dev.size(unit=c("cm"))
ggsave("Morphochisquare.svg", plot=Morphochisquarearranged, units="cm", width=50.80000, height=26.51125)
ggsave("Morphochisquare.tiff", plot=BladeletElwadarranged, units="cm", width=50.80000, height=26.51125, dpi=600, compression = 'lzw')