###### Community composition analysis #####
#
#Year: 2024

#The code runs a PERMANOVA to test differences in community composition between areas with and without _Lanice conchilega_
#It also runs a constrained analysis of proximities to visualize these differences, and a SIMPER analysis to identify species  
#which species contributed to these dissimilarities.

# R version used: 4.2.0
# Rstudio version used: RStudio 2023.06.2+561 "Mountain Hydrangea" Release (de44a3118f7963972e24a78b7a1ad48b4be8a217, 2023-08-23) for windows Mozilla/5.0 (Windows NT 10.0; Win64; x64) AppleWebKit/537.36 (KHTML, like Gecko) RStudio/2023.06.2+561 Chrome/110.0.5481.208 Electron/23.3.0 Safari/537.36

library(dplyr) # version 1.1.2
library(tidyr) # version 1.3.0
library(tibble) # version 3.2.1
library(vegan) # version 2.6-2
library(ggplot2) # version 3.4.2


##### Data ####
spp <- read.csv("code/community_composition_data.csv")


#### Permanova ####

#drop L. conchilega from the rest of the community matrix
spp1 <- spp %>% 
  dplyr::select(-sp_42)

#presence absence of L. conchilega
lc <- spp %>% 
  dplyr::select(sp_42) %>% 
  mutate(sp_42 = as.factor(sp_42))


##### Check variance ####
#calculate distance matrix
abun_bray <- vegdist(spp1, method = "bray")
#calculate variance
abun_bray_variance <- betadisper(abun_bray, lc$sp_42)

boxplot(abun_bray_variance)
anova(abun_bray_variance)#it is significant which might diminish the power of the test.

##### PERMANOVA ####
abun_bray_adonis <- adonis2(abun_bray ~ lc$sp_42)
abun_bray_adonis

##### SIMPER ####
sim <- simper(spp1, group = lc$sp_42)

##### Constrained analysis of proximities ####
constrained <- capscale(spp1 ~ lc$sp_42, distance = "bray")

contrained.scores <- as.data.frame(scores(constrained)$sites) %>% 
  mutate(sp_42 = lc$sp_42)

ggplot(contrained.scores, aes(CAP1,MDS1, color = sp_42))+
  geom_point()