ks_05_timeseries_mca-hcpc.Rmd

---
title: "ks_mca-hcpc_timeseries"
author: "Kate Sheridan"
date: "3/13/2022"
output: html_document
editor_options: 
  markdown: 
    wrap: 72
---

```{r setup, include=FALSE}
library(tidyverse)
library(ggplot2)
library(FactoMineR)
library(factoextra)
library(dplyr)
library(here)

# here filepaths
traits <- here('data', 'clean')
plots <- here('output', '2022may', 'mca')

# set seed for whole session:
addTaskCallback(function(...) {set.seed(3824);TRUE})
# to remove
#removeTaskCallback(1)
```

## load in data

```{r}
# time series
algae_traits <- read_csv(here(traits, "20220525_1982-end_algae.csv")) %>%
  filter(intertidal_gom == 'yes') %>%
  select(!(c(intertidal_gom, group, autotroph, 
             saccate_external_morphology,
             potential_asexual_reproduction,
             epilithic_env_position)))

#algae_traits <- algae_traits %>%
#  mutate(maximun_longevity = str_replace_all(maximun_longevity, c(" \\(one-ten years\\)"= "")))
  

animal_traits <- read_csv(here(traits, "20220321_1982-end_animal.csv")) %>%
  filter(intertidal == 'yes') %>%
  select(!(c(intertidal, trophic_level, pelagic, deposit_feeder, benthic)))

# to load in individual species lists from parts of the time series
#animal_1982_1995_traits <- read_csv(here(mcadata, "20220313_1982-1995_animal.csv")) %>%
#  filter(intertidal == 'yes') %>%
#  select(!(c(intertidal)))
#animal_2011_end_traits <- read_csv(here(mcadata, "20220313_2011-end_animal.csv")) %>%
#  filter(intertidal == 'yes') %>%
#  select(!(c(intertidal)))


# all species
#algae_imputed_all <- read_csv(here(traits, "20220323_algae_imputed.csv"),
#                                 col_select = !1) %>%
#  filter(intertidal == 'yes') %>%
#  select(!(c(intertidal)))

#animal_imputed_all <- read_csv(here(traits, "20220321_animal_imputed.csv"),
#                                 col_select = !1) %>%
#  filter(intertidal == 'yes') %>%
#  select(!(c(intertidal)))
```

## load in functions

```{r}
# input should be imputed traits, filtered as needed
## note that we're selecting columns automatically here 
## starting with column 3 after the species column; change if not true
# note that only output starting with gg are ggplot objects
appledoreMCAHCPC <- function(traits) {
  # select columns for MCA
  traits <- traits %>% 
    column_to_rownames('species') %>%
    dplyr::select(3:(ncol(traits)-1)) %>%
    mutate(across(.fns = ~ as.factor(.)))
  # run MCA
  mca2 <- MCA(traits, graph = FALSE)
  # extract variable categories
  mca2_vars <- get_mca_var(mca2)
  # prepare data for plotting
  cats=apply(traits, 2, function(x) nlevels(as.factor(x)))
  mca2_vars_df = data.frame(mca2$var$coord, Variable = rep(names(cats), 
                                                         cats))
  mca2_obs_df = data.frame(mca2$ind$coord)
  # plot with lines and vectors
  mcaggplot <- ggplot(data = mca2_obs_df, aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray70") + 
    geom_vline(xintercept = 0, colour = "gray70") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray80") + 
    geom_text(data = mca2_vars_df, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(mca2_vars_df), 
                                     colour = Variable)) + 
    ggtitle("MCA plot with FactorMineR") + 
    scale_colour_discrete(name = "Variable") +
    geom_segment(data = mca2_vars_df, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
    theme_bw(20)
  # ellipses plot
  plot_ellipses <- plotellipses(mca2)
  # diagnostics
  scree2 <- fviz_screeplot(mca2)
  contrib2 <- fviz_contrib(mca2, choice = "var", axes = 1:2, top = 15)
  cos2plot <- fviz_cos2(mca2, choice = 'var', axes = 1:2, top = 15)
  # this plot combines output of contribution and cos2, 
  ## eg: orange+high alpha = high contribution, high quality
  cos2contribplot <- fviz_mca_var(mca2, col.var = 'cos2', alpha.var = 'contrib',
                           gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"), 
                           repel = TRUE, # Avoid text overlapping
                           ggtheme = theme_minimal())
  
  
  # run HCPC
  hcpc2 <- HCPC(mca2, graph = FALSE)
  # hcpc dendrogram
  hcpc_dend <- fviz_dend(hcpc2, cex = .6)
  # hcpc factor map
  hcpc_clust <- fviz_cluster(hcpc2, geom = "point", main = "Factor map")
  
  
  #output
  list(mcares = mca2, mcacos2 = mca2_vars$cos2,
       vars_df = mca2_vars_df,
       clustmembers = hcpc2$desc.ind$para,
       clustervariables = hcpc2$desc.var$category,
       hcpc_chi = hcpc2$desc.var$test.chi2,
       ggmca = mcaggplot,
       ggscreeplot = scree2, ggcontrib = contrib2,
       ggcos2 = cos2plot, ggvars = cos2contribplot,
       ellipses = plot_ellipses, ggdendrogram = hcpc_dend,
       ggclusterplot = hcpc_clust)
}
```

## run function example

explore output; use df\$ to see all saved plots and summary stats Any
output starting with gg is a ggplot object and can be chained with
ggplot functions such as: `algae_all$ggmca + theme_bw()`

```{r function-explore}
algae_all_test <- algae_traits %>%
  select(!(c(year, transect, position))) %>%
  distinct()

algae_all <- appledoreMCAHCPC(algae_all_test)
# example view a summary
algae_all$mcares$eig
# example ggplot
algae_all$ggscreeplot + theme_bw()

algae_all$ggmca


algae_all$ggdendrogram
algae_all$ggclusterplot

algae_all$clustervariables
algae_all$clustmembers

```

# set up subsets and run

## Algae

```{r algae-all}
algaehcpc <- HCPC(algae_all$mcares, nb.clust = 5)

algaehcpc$desc.var
algaehcpc$desc.ind

fviz_dend(algaehcpc)

```

### Timeseries

```{r algae-timeseries}
# filter
algae_82_95 <- algae_traits %>%
  filter(year <= 1994 & year >= 1982) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

algae_96_06 <- algae_traits %>%
  filter(year >= 1996 & year < 2007) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

algae_10_end <- algae_traits %>%
  filter(year <= 2011) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

# run mca
mca_algae_82_95 <- appledoreMCAHCPC(algae_82_95)
mca_algae_96_06 <- appledoreMCAHCPC(algae_96_06)
mca_algae_10_end <- appledoreMCAHCPC(algae_10_end)
```

```{r algae-timeseries-explore}
mca_algae_82_95$ggmca
mca_algae_96_06$ggmca
mca_algae_10_end$ggmca
```

### exposed/sheltered

for the time series here it will just be beginning and end

```{r algae-exposedsheltered}
# just exposed v sheltered
algae_exposed <- algae_traits %>%
  filter(position == 'exposed') %>%
  select(!(c(year, transect, position))) %>%
  distinct()

algae_sheltered <- algae_traits %>%
  filter(position == 'sheltered') %>%
  select(!(c(year, transect, position))) %>%
  distinct()

mca_algae_exposed <- appledoreMCAHCPC(algae_exposed)
mca_algae_sheltered <- appledoreMCAHCPC(algae_sheltered)

# exposed v sheltered beginning and end
algae_exposed_start <- algae_traits %>%
  filter(position == 'exposed' & year <= 1994 & year >= 1982) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

algae_sheltered_start <- algae_traits %>%
  filter(position == 'sheltered' & year <= 1994 & year >= 1982) %>%
  select(!(c(year, transect, position))) %>%
  distinct()


algae_sheltered_middle <- algae_traits %>%
  filter(position == 'sheltered' & year >= 1996 & year < 2007) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

algae_exposed_middle <- algae_traits %>%
  filter(position == 'exposed' & year >= 1996 & year < 2007) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

algae_exposed_end <- algae_traits %>%
  filter(position == 'exposed' & year <= 2010) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

algae_sheltered_end <- algae_traits %>%
  filter(position == 'sheltered' & year <= 2010) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

mca_algae_exposed_start <- appledoreMCAHCPC(algae_exposed_start)
mca_algae_sheltered_start <- appledoreMCAHCPC(algae_sheltered_start)
mca_algae_exposed_middle <- appledoreMCAHCPC(algae_exposed_middle)
mca_algae_sheltered_middle <- appledoreMCAHCPC(algae_sheltered_middle)
mca_algae_exposed_end <- appledoreMCAHCPC(algae_exposed_end)
mca_algae_sheltered_end <- appledoreMCAHCPC(algae_sheltered_end)
```

```{r algae-position-plots}

mca_algae_exposed_start$ggmca
mca_algae_exposed_end$ggmca
mca_algae_sheltered_start$ggmca
mca_algae_sheltered_end$ggmca
mca_algae_exposed$ggmca
mca_algae_sheltered$ggmca

ggsave(here(plots, 'temporary', '2algaeexposed.png'), mca_algae_exposed$ggmca, width = 10)
ggsave(here(plots, 'temporary', '2algaesheltered.png'), mca_algae_sheltered$ggmca, width = 10)

ggsave(here(plots, 'temporary', '2algaeexposed_start.png'), mca_algae_exposed_start$ggmca, width = 10)
ggsave(here(plots, 'temporary', '2algaeexposed_end.png'), mca_algae_exposed_end$ggmca, width = 10)

ggsave(here(plots, 'temporary', '2algaesheltered_start.png'), mca_algae_sheltered_start$ggmca, width = 10)
ggsave(here(plots, 'temporary', '2algaesheltered_end.png'), mca_algae_sheltered_end$ggmca, width = 10)

```

```{r algae-exposedsheltered-explore}


mca_algae_exposed$ggdendrogram
mca_algae_sheltered$ggdendrogram

mca_algae_exposed_start$ggdendrogram
mca_algae_exposed_end$ggdendrogram

mca_algae_sheltered_start$ggdendrogram
mca_algae_sheltered_end$ggdendrogram

mca_algae_sheltered_end$clustmembers
mca_algae_sheltered_end$clustervariables

mca_algae_sheltered_start$clustmembers
mca_algae_sheltered_start$clustervariables


exposedhcpc2 <- HCPC(res = mca_algae_exposed_end$mcares, nb.clust =  5)
mca_algae_exposed_start$ggdendrogram
fviz_dend(exposedhcpc2)







```

```{r algae-exposedsheltered-summarystats}
mca_algae_exposed_start$mcares$eig
mca_algae_exposed_end$mcares$eig
mca_algae_exposed_middle$mcares$eig
mca_algae_sheltered_start$mcares$eig
mca_algae_exposed_middle$mcares$eig
mca_algae_sheltered_end$mcares$eig


mca_algae_exposed_start$ggcos2
mca_algae_exposed_middle$ggcos2
mca_algae_exposed_end$ggcos2
mca_algae_sheltered_start$ggcos2
mca_algae_sheltered_middle$ggcos2
mca_algae_sheltered_end$ggcos2


mca_algae_exposed_start$mcares$var$contrib


mca_algae_exposed_start$ggcontrib
mca_algae_exposed_middle$ggcontrib
mca_algae_exposed_end$ggcontrib
mca_algae_sheltered_start$ggcontrib
mca_algae_sheltered_middle$ggcontrib
mca_algae_sheltered_end$ggcontrib


mca_algae_sheltered_end$mcares$ind$contrib
mca_algae_sheltered_end$mcares$ind$cos2

mca_algae_sheltered_end$mcares$var$eta2
mca_algae_sheltered_end$mcares$var$contrib

mca_algae_sheltered_start$mcares$ind$contrib
mca_algae_sheltered_start$mcares$ind$cos2

mca_algae_sheltered_start$mcares$var$eta2
mca_algae_sheltered_start$mcares$var$contrib

```

```{r algae-compare}
library(arsenal)
# in output, first dataframe is X, second is Y
exposedvprotect <- comparedf(algae_exposed, algae_sheltered, by = 'species')
#extract differences
diffs(exposedvprotect, what = c('observations'))

protectstvend <- comparedf(algae_sheltered_start, algae_sheltered_end, by = 'species')
diffs(protectstvend, what = c('observations'))

exposedstvend <- comparedf(algae_exposed_start, algae_exposed_end, by = 'species')
diffs(exposedstvend, what = c('observations'))

```


```{r algae-hypervolume}
library(hypervolume)
#testing
#hypervolume1<-hypervolume_box(data.frame(algae_all$mcares$ind$coord))
#plot(hypervolume1)


hyp_algae_exposed_start <-hypervolume_gaussian(data.frame(mca_algae_exposed_start$mcares$ind$coord))
#hyp_algae_exposed_end <-hypervolume_box(data.frame(mca_algae_exposed_end$mcares$ind$coord))
hyp_algae_sheltered_start <-hypervolume_gaussian(data.frame(mca_algae_sheltered_start$mcares$ind$coord))
hyp_algae_sheltered_end <-hypervolume_gaussian(data.frame(mca_algae_sheltered_end$mcares$ind$coord))


# protected v exposed start
hyp_set_algae_sides1 <- hypervolume_set(hyp_algae_exposed_start, hyp_algae_sheltered_start, check.memory=FALSE)
# 1 is first listed, 2 is second listed
hypervolume_overlap_statistics(hyp_set_algae_sides1)
#volumes
get_volume(hyp_set_algae_sides1)
# basic plot
plot.HypervolumeList(hyp_set_algae_sides1, colors = c('#806238', "#94b4d8", '#bcbcd2', 'purple', '#a08b5a', 'blue'))


# protected v exposed current
hyp_set_algae_sides <- hypervolume_set(hyp_algae_exposed_end, hyp_algae_sheltered_end, check.memory=FALSE)
# 1 is first listed, 2 is second listed
hypervolume_overlap_statistics(hyp_set_algae_sides)
#volumes
get_volume(hyp_set_algae_sides)
# basic plot
plot.HypervolumeList(hyp_set_algae_sides, colors = c('#806238', "#94b4d8", '#bcbcd2', 'purple', '#a08b5a', 'blue'))



# exposed beginning v end
hyp_set_algae_exposed <- hypervolume_set(hyp_algae_exposed_start, hyp_algae_exposed_end, check.memory=FALSE)
# 1 is first listed, 2 is second listed
hypervolume_overlap_statistics(hyp_set_algae_exposed)
#volumes
get_volume(hyp_set_algae_exposed)
# basic plot
plot.HypervolumeList(hyp_set_algae_exposed, colors = c('#806238', "#94b4d8", '#bcbcd2', 'purple', '#a08b5a', 'blue'))


# protected beginning v end
hyp_set_algae_sheltered <- hypervolume_set(hyp_algae_sheltered_start, hyp_algae_sheltered_end, check.memory=FALSE)
# 1 is first listed, 2 is second listed
hypervolume_overlap_statistics(hyp_set_algae_sheltered)
#volumes
get_volume(hyp_set_algae_sheltered)
# basic plot
plot.HypervolumeList(hyp_set_algae_sheltered, colors = c('#806238', "#94b4d8", '#bcbcd2', 'purple', '#a08b5a', 'blue'),
                     show.axes = FALSE, cex.random = .1)

hypervolume_redundancy(hyp_algae_sheltered_start, data.frame(mca_algae_sheltered_end$mcares$ind$coord))
hypervolume_redundancy(hyp_algae_sheltered_start, data.frame(mca_algae_sheltered_start$mcares$ind$coord))


## both tests seem to show change, but not significance. 

# to run this test I need to sample from both populations not just one; 
# run resample on all sheltered, then overlap test with the path being to overlap test
## this method is less good for small sample sizes
hyp_algae_shelter_all <- hypervolume_gaussian(data.frame(mca_algae_sheltered_end$mcares$ind$coord))
test_algae_start <- hypervolume_resample('test_algae_start', hyp_algae_shelter_all, 'bootstrap', cores = 2)

# this is good for smaller datasets
test_algae_start <- hypervolume_permute('test_algae_permute', hyp_algae_sheltered_start, hyp_algae_sheltered_end, cores = 2)
algae_sheltered_overlap <- hypervolume_overlap_test(hyp_algae_sheltered_start, hyp_algae_sheltered_end, path = test_algae_start)

algae_sheltered_overlap

hypervolume_funnel()

plot(hyp_algae_sheltered_start)
plot(hyp_algae_sheltered_end)
```




## Animals

Note animal_all gives different result than any group of animals. I
really don't understand why, it puts Astarias rubens in totally
different locations.

```{r animal-all}
animal_all_t <- animal_traits %>%
  select(!(c(year,transect,position))) %>%
  distinct()

mca_animal_all <- appledoreMCAHCPC(animal_all_t)

animal_allhcpc <- HCPC(res = mca_animal_all$mcares, nb.clust =  5)

# plots testing output
fviz_cluster(animal_allhcpc, geom = "point", main = "Factor map")
mca_animal_all$ggclusterplot
fviz_dend(animal_allhcpc)
mca_animal_all$clustmembers
mca_animal_all$clustervariables

animal_allhcpc$desc.var
animal_allhcpc$desc.ind
```

### Timeseries

```{r animal-timeseries}
# filter
animal_82_95 <- animal_traits %>%
  filter(year <= 1994 & year >= 1982) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

animal_96_06 <- animal_traits %>%
  filter(year >= 1996 & year < 2007) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

animal_10_end <- animal_traits %>%
  filter(year <= 2011) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

# run mca
mca_animal_82_95 <- appledoreMCAHCPC(animal_82_95)
mca_animal_96_06_mca <- appledoreMCAHCPC(animal_96_06)
mca_animal_10_end <- appledoreMCAHCPC(animal_10_end)
```

### exposed/sheltered

for the time series here it will just be beginning and end

```{r animal-exposedsheltered}
# just exposed v sheltered
animal_exposed <- animal_traits %>%
  filter(position == 'exposed') %>%
  select(!(c(year, transect, position))) %>%
  distinct()

animal_sheltered <- animal_traits %>%
  filter(position == 'sheltered') %>%
  select(!(c(year, transect, position))) %>%
  distinct()

mca_animal_exposed <- appledoreMCAHCPC(animal_exposed)
mca_animal_sheltered <- appledoreMCAHCPC(animal_sheltered)

# exposed v sheltered beginning and end
animal_exposed_start <- animal_traits %>%
  filter(position == 'exposed' & year <= 1994 & year >= 1982) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

animal_sheltered_start <- animal_traits %>%
  filter(position == 'sheltered' & year <= 1994 & year >= 1982) %>%
  select(!(c(year, transect, position))) %>%
  distinct()


animal_exposed_middle <- animal_traits %>%
  filter(position == 'exposed' & year >= 1996 & year < 2007) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

animal_sheltered_middle <- animal_traits %>%
  filter(position == 'sheltered' & year >= 1996 & year < 2007) %>%
  select(!(c(year, transect, position))) %>%
  distinct()


animal_exposed_end <- animal_traits %>%
  filter(position == 'exposed' & year <= 2010) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

animal_sheltered_end <- animal_traits %>%
  filter(position == 'sheltered' & year <= 2010) %>%
  select(!(c(year, transect, position))) %>%
  distinct()

mca_animal_exposed_start <- appledoreMCAHCPC(animal_exposed_start)
mca_animal_sheltered_start <- appledoreMCAHCPC(animal_sheltered_start)
mca_animal_exposed_middle <- appledoreMCAHCPC(animal_exposed_middle)
mca_animal_sheltered_middle <- appledoreMCAHCPC(animal_sheltered_middle)
mca_animal_exposed_end <- appledoreMCAHCPC(animal_exposed_end)
mca_animal_sheltered_end <- appledoreMCAHCPC(animal_sheltered_end)
```

```{r animal-position-plots}
# overall
mca_animal_exposed$ggmca
mca_animal_sheltered$ggmca


mca_animal_exposed_start$ggvars
mca_animal_exposed_end$ggvars


# exposed
mca_animal_exposed_start$ggmca
mca_animal_exposed_end$ggmca


# sheltered
mca_animal_sheltered_start$ggmca
mca_animal_sheltered_end$ggmca

mca_animal_sheltered_start$ggvars
mca_animal_sheltered_end$ggvars

ggsave(here(plots, 'temporary', 'animalexposed_start.png'), mca_animal_exposed_start$ggmca, width = 10)
ggsave(here(plots, 'temporary', 'animalexposed_end.png'), mca_animal_exposed_end$ggmca, width = 10)

ggsave(here(plots, 'temporary', 'animalsheltered_start.png'), mca_animal_sheltered_start$ggmca, width = 10)
ggsave(here(plots, 'temporary', 'animalsheltered_end.png'), mca_animal_sheltered_end$ggmca, width = 10)

```

```{r animal-hcpc-exploration}
shelteredhcpc2 <- HCPC(res = mca_animal_sheltered$mcares, nb.clust =  5)
exposedhcpc2 <- HCPC(res = mca_animal_exposed$mcares, nb.clust =  4)
  
fviz_cluster(exposedhcpc2, geom = "point", main = "Factor map") 
fviz_cluster(shelteredhcpc2, geom = "point", main = "Factor map")

fviz_dend(exposedhcpc2)
fviz_dend(shelteredhcpc2)

mca_animal_exposed$ggclusterplot
mca_animal_sheltered$ggclusterplot

mca_animal_exposed$ggdendrogram
mca_animal_sheltered$ggdendrogram

mca_animal_exposed$clustervariables
mca_animal_sheltered$clustervariables

#mca_animal_exposed$clustmembers
mca_animal_sheltered$clustmembers

mca_animal_sheltered_start$ggclusterplot
mca_animal_sheltered_end$ggclusterplot
mca_animal_exposed_start$ggclusterplot
mca_animal_exposed_end$ggclusterplot

mca_animal_sheltered$ggdendrogram
mca_animal_sheltered_start$ggdendrogram
mca_animal_sheltered_end$ggdendrogram

mca_animal_exposed$ggdendrogram
mca_animal_exposed_start$ggdendrogram
mca_animal_exposed_end$ggdendrogram


mca_animal_sheltered_end$mcares$ind$cos2
mca_animal_sheltered_end$mcares$ind$contrib

mca_animal_sheltered_end$ggvars


fviz_mca_ind(mca_animal_sheltered_end$mcares, col.var = 'contrib', #alpha.var = 'contrib',
                           gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
                           select.var = list(cos2 = .4),
                           repel = TRUE, # Avoid text overlapping
                           ggtheme = theme_minimal()) 

mca_animal_exposed_end$mcares$var$eta2
mca_animal_exposed_end$mcares$var$contrib



```

# arsenal to extract differences quickly

```{r animal-compare}
library(arsenal)
# in output, first dataframe is X, second is Y
exposedvprotect <- comparedf(animal_exposed, animal_sheltered, by = 'species')
#extract differences
diffs(exposedvprotect, what = c('observations'))

protectstvend <- comparedf(animal_sheltered_start, animal_sheltered_end, by = 'species')
diffs(protectstvend, what = c('observations'))

exposedstvend <- comparedf(animal_exposed_start, animal_exposed_end, by = 'species')
diffs(exposedstvend, what = c('observations'))


animal_exp_notmiddle <- full_join(animal_exposed_start, animal_exposed_end)

exposedstvend <- comparedf(animal_exposed_middle, animal_exp_notmiddle, by = 'species')
diffs(exposedstvend, what = c('observations'))

#nothing unique to the middle

```


```{r animal-hypervolume}
library(hypervolume)

hyp_animal_sheltered_start <-hypervolume_gaussian(data.frame(mca_animal_sheltered_start$mcares$ind$coord))
hyp_animal_sheltered_end <-hypervolume_gaussian(data.frame(mca_animal_sheltered_end$mcares$ind$coord))

# protected beginning v end
hyp_set_animal_sheltered <- hypervolume_set(hyp_animal_sheltered_start, hyp_animal_sheltered_end, check.memory=FALSE)
# 1 is first listed, 2 is second listed
hypervolume_overlap_statistics(hyp_set_animal_sheltered)
#volumes
get_volume(hyp_set_animal_sheltered)
# basic plot
plot.HypervolumeList(hyp_set_animal_sheltered, colors = c('#806238', "#94b4d8", '#bcbcd2', 'purple', '#a08b5a', 'blue'),
                     show.axes = FALSE, cex.random = .1)


```

# presentation plots

palette is:
`#bcbcd2, #94b4d8, #5484c4, #374658, #a08b5a, #806238, #543414, #3a230e`
lavender grey, pale cerulean, glaucous, charcoal, metallic sunburst,
coyote brown, dark brown, bistre

sets of 2 should be primarily drawn from the middle two of each part of
the palette: 
two blues 
#94b4d8, #5484c4 
two browns 
#806238, #543414 
blue brown light 
#94b4d8,#806238 
blue brown dark 
#5484c4, #543414

```{r nice-fonts}
library(showtext)
showtext_auto()
# to see list of fonts
font_files()
# add desired fonts with ('name to call it by in code', 'name in filepath, or complete filepath')
font_add('quicksand', 'Quicksand-Medium.ttf')
```



## clusters

```{r testing}
#tree
library(dendextend)

# extract dendrogram and save as ggdend object
test <- as.ggdend(mca_animal_exposed$ggdendrogram$plot_env$dend)

ggplot(test, horiz = TRUE) + 
  labs(y = NULL, title = 'Algae Clusters') +
  ylim(.25,-.12) +
  # manually set colors
  scale_color_manual(values = c('red','blue','purple','black','green','orange','grey')) +
  theme(axis.ticks = element_blank(),
        axis.text.y = element_blank()
        )

ggsave(filename = here(plots, 'algaetest.png'), device = 'png', 
       height = 5, width = 7, units = 'in')



# clusters

mca_algae_exposed_start$ggclusterplot +
  #manually set colors
  scale_fill_manual(values = c('red','blue','purple','black','green','orange')) +
  scale_color_manual(values = c('red','blue','purple','black','green','orange')) +
  labs(title = 'Algae clusters') +
  theme_classic() +
  theme(axis.text = element_blank())

ggsave(filename = here(plots, 'algaetest2.png'), device = 'png', 
       height = 5, width = 7, units = 'in')


#mcas
mca_animal_exposed_start$ggmca + 
  labs(title = 'Algae mca') +
    scale_color_manual(values = c('red','blue','purple', 'black',
                                  'green','orange', 'grey', 'yellow', 'brown',
                                  'grey17', 'pink', 'grey60', 'salmon',
                                  'violet', 'cornsilk2')) +
  theme_classic() +
  theme(legend.position = 'none')
  
  



  
ggsave(filename = here(plots, 'algaetest3.png'), device = 'png', 
       height = 5, width = 7, units = 'in')

```

```{r algae}
library(dendextend)
algaehcpc_dend <- fviz_dend(algaehcpc, cex = .6)
algaehcpc_dend <- as.ggdend(algaehcpc_dend$plot_env$dend)
algaehcpc_cluster <- fviz_cluster(algaehcpc, geom = "point", main = "Factor map")


ggplot(algaehcpc_dend, horiz = TRUE) + 
  labs(y = NULL, title = 'Algae Functional Groups: All') +
  ylim(.2,-.12) +
  # manually set colors
  scale_color_manual(values = c('#bcbcd2','#5484c4','#543414','#374658', '#a08b5a', 'black')) +
  theme(axis.ticks = element_blank(),
        axis.text.y = element_blank()
        )

ggsave(filename = here(plots, 'hcpc', '20220525_algae_dend.png'), device = 'png', 
       height = 5, width = 7, units = 'in')

algaehcpc_cluster$data


algaehcpc_cluster +
  #manually set colors
  scale_fill_manual(values = c('#543414', '#a08b5a','#5484c4', '#bcbcd2','#374658')) +
  scale_color_manual(values = c('#543414','#a08b5a','#5484c4', '#bcbcd2','#374658')) +
  labs(title = 'Algae functional groups: All') +
  theme_classic() +
  theme(axis.text = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, 'hcpc', '20220323_algae_all_cluster.png'), device = 'png', 
       height = 5, width = 7, units = 'in')

```



```{r animal-all}
animal_allhcpc <- HCPC(res = mca_animal_all$mcares, nb.clust =  4)
animal_hcpc_dend <- fviz_dend(animal_allhcpc, cex = .6)
animal_hcpc_dend <- as.ggdend(animal_hcpc_dend$plot_env$dend)
animalhcpc_cluster <- fviz_cluster(animal_allhcpc, geom = "point", main = "Factor map")

ggplot(animal_hcpc_dend, horiz = TRUE) + 
  labs(y = NULL, title = 'Animal Functional groups: All') +
  ylim(.34,-.12) +
  # manually set colors
  scale_color_manual(values = c('#5484c4','#543414','#374658', '#a08b5a', 'black')) +
  theme(axis.ticks = element_blank(),
        axis.text.y = element_blank()
        )

ggsave(filename = here(plots, 'hcpc', '20220323_animal_dend_all.png'), device = 'png', 
       height = 5, width = 7, units = 'in')

animalhcpc_cluster$data

animalhcpc_cluster +
  #manually set colors
  scale_fill_manual(values = c('#374658', '#5484c4', '#a08b5a', '#543414')) +
  scale_color_manual(values = c('#374658','#5484c4', '#a08b5a', '#543414')) +
  labs(title = 'Animal functional groups: All') +
  theme_classic() +
  theme(axis.text = element_blank(),
        legend.position = 'none')


ggsave(filename = here(plots, 'hcpc', '20220323_animal_all_cluster.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
```



```{r animals}
#manually split the 4 clusters for sheltered_start
## probably better plotting can be done at this point but 
# i'm just going to make it into the same format as the others
shelteredhcpc2 <- HCPC(res = mca_animal_sheltered_start$mcares, nb.clust =  5)
start_test <- fviz_dend(shelteredhcpc2, cex = .7)
start_test2 <- as.ggdend(start_test$plot_env$dend)
start_cluster2 <- fviz_cluster(shelteredhcpc2, geom = "point", main = "Factor map")


shelteredhcpc3 <- HCPC(res = mca_animal_sheltered_end$mcares, nb.clust =  5)
end_test <- fviz_dend(shelteredhcpc3, cex = .7)
end_test2 <- as.ggdend(end_test$plot_env$dend)


#dend_animal_sheltered_end <- as.ggdend(mca_animal_sheltered_end$ggdendrogram$plot_env$dend)


ggplot(start_test2, horiz = TRUE) + 
  labs(y = NULL, title = 'Animal Functional groups: Protected 1982-1995') +
  ylim(.25,-.12) +
  # manually set colors
  scale_color_manual(values = c('#5484c4','#bcbcd2','#543414','#374658', '#a08b5a', 'black')) +
  theme(axis.ticks = element_blank(),
        axis.text.y = element_blank()
        )

#ggsave(filename = here(plots, 'hcpc', '20220321_animal_dend_protected-start.png'), device = 'png', 
#       height = 5, width = 7, units = 'in')


ggplot(end_test2, horiz = TRUE) + 
  labs(y = NULL, title = 'Animal functional groups: protected 2010-present') +
  ylim(.3,-.12) +
  # manually set colors
  scale_color_manual(values = c('#a08b5a','#543414','#bcbcd2','#5484c4', '#374658', 'black')) +
  theme(axis.ticks = element_blank(),
        axis.text.y = element_blank()
        )

ggsave(filename = here(plots, 'hcpc', '20220321_animal_protected-end_dend.png'), device = 'png', 
       height = 5, width = 7, units = 'in')



start_cluster2 +
  #manually set colors
  scale_fill_manual(values = c('#543414', '#a08b5a', '#5484c4', '#bcbcd2')) +
  scale_color_manual(values = c('#543414', '#a08b5a', '#5484c4', '#bcbcd2')) +
  labs(title = 'Animal functional groups: protected 1982-1995') +
  theme_classic() +
  theme(axis.text = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, 'hcpc', '20220321_animal_protected-start_cluster.png'), device = 'png', 
       height = 5, width = 7, units = 'in')

mca_animal_sheltered_end$ggclusterplot +
  #manually set colors
  scale_fill_manual(values = c('#a08b5a','#543414','#bcbcd2')) +
  scale_color_manual(values = c('#a08b5a','#543414','#bcbcd2')) +
  labs(title = 'Animal functional groups: protected 2010-present') +
  theme_classic() +
  theme(axis.text = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, 'hcpc', '20220321_animal_protected-end_cluster.png'), device = 'png', 
       height = 5, width = 7, units = 'in')





dend_animal_exposed <- as.ggdend(mca_animal_exposed_end$ggdendrogram$plot_env$dend)

ggplot(dend_animal_exposed, horiz = TRUE) + 
  labs(y = NULL, title = 'Animal functional groups: exposed 2010-present') +
  ylim(.3,-.12) +
  # manually set colors
  scale_color_manual(values = c('#a08b5a','#543414','#bcbcd2','black')) +
  theme(axis.ticks = element_blank(),
        axis.text.y = element_blank()
        )

ggsave(filename = here(plots, 'hcpc', '20220321_animal_exposed-end_dend.png'), device = 'png', 
       height = 5, width = 7, units = 'in')

```

## MCA plots


```{r algae-scratch}


fviz_mca_var(mca_animal_sheltered_end$mcares, col.var = 'cos2', alpha.var = 'contrib',
                           gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
                           select.var = list(cos2 = .4),
                           repel = TRUE, # Avoid text overlapping
                           ggtheme = theme_minimal()) 


mca_algae_sheltered_middle$vars_df

#'space_use_sub_canopy', 'crustose_morphology_growth_form', 'leathery_external_morphology', 'erect_morphology_growth_form'

#('branched_external_morphology', 'leathery_external_morphology', 'erect_morphology_growth_form')
```

```{r algae-all-mca}
library(ggrepel)
#added scavenger and filter/suspension to >.6 cos2
# make subset of vars to plot arrows for
arrows1 <- algae_all$vars_df %>%
  filter(Variable %in% c('branched_external_morphology', 'leathery_external_morphology')) %>%
  `row.names<-`(c('Not Leathery', 'Leathery', 'Not Branched', 'Branched'))

# exposed
ggplot(data.frame(algae_all$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = arrows1, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(arrows1), 
                                     colour = Variable), seed = 3824) +
  labs(title = "Algae: all", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#5484c4','purple')) +
    geom_segment(data = arrows1, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.25,1.8) +
  ylim(-1.15,1.7) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220525_algae_all.png'), device = 'png', 
       height = 5, width = 7, units = 'in')


```

```{r}
ggplot() +
  geom_hline(yintercept = 0, colour = "gray90") + 
  geom_vline(xintercept = 0, colour = "gray90") + 
  geom_density2d(data = as.data.frame(mca_algae_sheltered_start$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "#94b4d8") +
  #geom_density2d(data = as.data.frame(mca_animal_sheltered_middle$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "grey60") +
  geom_density2d(data = as.data.frame(mca_algae_sheltered_end$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "#374658") +
  ylim(-1.2,1.3) +
  xlim(-1.5,1.5) +
  labs(title = 'Algae: Protected', subtitle = 'Beginning (blue) to end (black)', x = NA, y = NA) +
  theme_classic() +
  theme(axis.text = element_blank())

ggsave(filename = here(plots, '20220525_algae_sheltered_change.png'), device = 'png', 
       height = 5, width = 7, units = 'in')


ggplot() +
  geom_hline(yintercept = 0, colour = "gray90") + 
  geom_vline(xintercept = 0, colour = "gray90") + 
  geom_density2d(data = as.data.frame(mca_algae_exposed_start$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "#94b4d8") +
  #geom_density2d(data = as.data.frame(mca_animal_sheltered_middle$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "grey60") +
  geom_density2d(data = as.data.frame(mca_algae_exposed_end$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "#374658") +
  ylim(-1.2,1.3) +
  xlim(-1.5,1.5) +
  labs(title = 'Algae: Exposed', subtitle = 'Beginning (blue) to end (black)', x = NA, y = NA) +
  theme_classic() +
  theme(axis.text = element_blank())

ggsave(filename = here(plots, '20220525_algae_exposed_change.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
```


```{r algae}
library(ggrepel)
#added scavenger and filter/suspension to >.6 cos2
# make subset of vars to plot arrows for
exposed_arrows1 <- mca_algae_exposed_start$vars_df %>%
  filter(Variable %in% c('branched_external_morphology', 'leathery_external_morphology')) %>%
  `row.names<-`(c('Not Leathery', 'Leathery', 'Not Branched', 'Branched'))

exposed_arrows2 <- mca_algae_exposed_middle$vars_df %>%
  filter(Variable %in% c('branched_external_morphology', 'leathery_external_morphology')) %>%
  `row.names<-`(c('Not Leathery', 'Leathery', 'Not Branched', 'Branched'))

exposed_arrows3 <- mca_algae_exposed_end$vars_df %>%
  filter(Variable %in% c('branched_external_morphology', 'leathery_external_morphology')) %>%
  `row.names<-`(c('Not Leathery', 'Leathery', 'Not Branched', 'Branched'))


# exposed
ggplot(data.frame(mca_algae_exposed_start$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = exposed_arrows1, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(exposed_arrows1), 
                                     colour = Variable,
                                     family = 'quicksand'
                                     ), seed = 3824) +
  labs(title = "Algae: 1982-1995, Exposed", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#5484c4','purple')) +
    geom_segment(data = exposed_arrows1, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.25,1.75) +
  ylim(-1.15,1.45) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220525_algae_exposed_start.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
#ggsave(filename = here(plots, '20220327_algae_exposed_start.pdf'), device = 'pdf', 
#       height = 5, width = 7, units = 'in')


ggplot(data.frame(mca_algae_exposed_middle$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = exposed_arrows2, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(exposed_arrows2), 
                                     colour = Variable,
                                     family = 'quicksand'
                                     ), seed = 3824) +
  labs(title = "Algae: 1996-2006, Exposed", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#5484c4','purple')) +
    geom_segment(data = exposed_arrows2, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.25,1.75) +
  ylim(-1.15,1.45) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220525_algae_exposed_mid.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
#ggsave(filename = here(plots, '20220327_algae_exposed_mid.pdf'), device = 'pdf', 
#       height = 5, width = 7, units = 'in')


ggplot(data.frame(mca_algae_exposed_end$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = exposed_arrows3, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(exposed_arrows3), 
                                     colour = Variable,
                                     family = 'quicksand'
                                     ), seed = 3824) +
  labs(title = "Algaes: 2010-present, Exposed", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#5484c4','purple')) +
    geom_segment(data = exposed_arrows3, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.25,1.75) +
  ylim(-1.15,1.45) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220525_algae_exposed_end.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
#ggsave(filename = here(plots, '20220327_algae_exposed_end.pdf'), device = 'pdf', 
#       height = 5, width = 7, units = 'in')




# make subset of vars to plot arrows for
sheltered_arrows2 <- mca_algae_sheltered_end$vars_df %>%
  filter(Variable %in% c('branched_external_morphology', 'filamentous_morphology_growth_form', 
                         'leathery_external_morphology', 'erect_morphology_growth_form')) %>%
  `row.names<-`(c('Not erect growth','Erect growth','Not Filamentous','Filamentous', 'Not Leathery', 'Leathery', 'Not Branched', 'Branched')) %>%
  filter(!(Variable == 'leathery_external_morphology' & Dim.1 < 0))



# sheltered
sheltered_arrows1 <- mca_algae_sheltered_start$vars_df %>%
  filter(Variable %in% c('branched_external_morphology', 'filamentous_morphology_growth_form', 
                         'leathery_external_morphology', 'erect_morphology_growth_form')) %>%
  `row.names<-`(c('Not erect growth','Erect growth','Not Filamentous','Filamentous', 'Not Leathery', 'Leathery', 'Not Branched', 'Branched')) #%>%
  #filter(!(Variable == 'scavenger' & Dim.1 < 0 | Variable == 'epibiotic' & Dim.1 > 0 | Variable == 'filter_suspension' & Dim.1 > 0))

#'space_use_sub_canopy', 'crustose_morphology_growth_form', 'leathery_external_morphology', 'erect_morphology_growth_form'

sheltered_arrows2 <- mca_algae_sheltered_middle$vars_df %>%
  filter(Variable %in% c('branched_external_morphology', 'crustose_morphology_growth_form', 
                         'leathery_external_morphology', 'erect_morphology_growth_form',
                         'space_use_sub_canopy')) %>%
  `row.names<-`(c('Not erect growth','Erect growth','Not Crustose','Crustose',
                  'Not Leathery', 'Leathery', 'Not Branched', 'Branched', 'Not Subcanopy', 'Subcanopy')) %>%
  filter(!(Variable == 'crustose_morphology_growth_form' & Dim.1 > 0 | Variable == 'erect_morphology_growth_form' & Dim.1 < 0))

sheltered_arrows3 <- mca_algae_sheltered_end$vars_df %>%
  filter(Variable %in% c('branched_external_morphology', 'filamentous_morphology_growth_form', 
                         'leathery_external_morphology', 'erect_morphology_growth_form')) %>%
  `row.names<-`(c('Not erect growth','Erect growth','Not Filamentous','Filamentous', 'Not Leathery', 'Leathery', 'Not Branched', 'Branched')) %>%
  filter(!(Variable == 'filamentous_morphology_growth_form' & Dim.1 > 0 | Variable == 'leathery_external_morphology' & Dim.1 < 0))


ggplot(data.frame(mca_algae_sheltered_start$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = sheltered_arrows1, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(sheltered_arrows1), 
                                     colour = Variable,
                                     family = 'quicksand'
                                     ), seed = 3824) +
  labs(title = "algaes: 1982-1994, Protected", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#374658', '#94b4d8', '#5484c4','orange', 'blue')) +
    geom_segment(data = sheltered_arrows1, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.25,1.75) +
  ylim(-1.15,1.45) +
  #scale_y_reverse(limits = c(1.3,-2.5)) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')


ggsave(filename = here(plots, '20220327_algae_sheltered-start.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
ggsave(filename = here(plots, '20220327_algae_sheltered-start.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')


ggplot(data.frame(mca_algae_sheltered_middle$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = sheltered_arrows2, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(sheltered_arrows2), 
                                     colour = Variable,
                                     family = 'quicksand'
                                     ), seed = 3824) +
  labs(title = "Algaes: 1996-2006, Protected", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#a08b5a', '#374658','#5484c4', '#543414', 'blue')) +
    geom_segment(data = sheltered_arrows2, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.25,1.75) +
  ylim(-1.15,1.45) +
  #scale_y_reverse(limits = c(1.3,-2.5)) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')


  
ggsave(filename = here(plots, '20220327_algae_sheltered-middle.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
ggsave(filename = here(plots, '20220327_algae_sheltered-middle.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')


ggplot(data.frame(mca_algae_sheltered_end$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = sheltered_arrows3, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(sheltered_arrows3), 
                                     colour = Variable,
                                     family = 'quicksand'
                                     ), seed = 3824) +
  labs(title = "Algae: 2010-present, Protected", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238', '#374658', '#94b4d8','#5484c4', 'orange', 'blue')) +
    geom_segment(data = sheltered_arrows3, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.25,1.75) +
  ylim(-1.15,1.45) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220327_algae_sheltered-end.png'), device = 'png', 
       height = 5, width = 7, units = 'in')  
ggsave(filename = here(plots, '20220327_algae_sheltered-end.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')  
  







```






```{r}
mca_animal_exposed_start$ggmca
mca_animal_exposed_end$ggmca


# sheltered
mca_animal_sheltered_start$ggmca
mca_animal_sheltered_end$ggmca

mca_animal_sheltered_start$ggvars
mca_animal_sheltered_end$ggvars


mca_animal_exposed_end$ggmca + 
  labs(title = 'Animals: Exposed') +
    scale_color_manual(values = c('red','blue','purple', 'black',
                                  'green','orange', 'grey', 'yellow', 'brown',
                                  'grey17', 'pink', 'grey60', 'salmon',
                                  'violet', 'cornsilk2')) +
  theme_classic() +
  theme(legend.position = 'none')

#how to make ggvars
fviz_mca_var(mca_animal_sheltered_end$mcares, col.var = 'cos2', alpha.var = 'contrib',
                           gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
                           select.var = list(cos2 = .6),
                           repel = TRUE, # Avoid text overlapping
                           ggtheme = theme_minimal()) +
  geom_density2d(colour = "gray80")


mca_animal_exposed_end$ggvars +
  geom_segment(data = mca_animal_exposed_end$vars_df, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  geom_density2d(colour = "gray80")

# how to mca plot
ggplot(data.frame(mca_animal_exposed_end$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray70") + 
    geom_vline(xintercept = 0, colour = "gray70") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray80") + 
    geom_text(data = mca_animal_exposed_end$vars_df, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(mca_animal_exposed_end$vars_df), 
                                     colour = Variable)) + 
    ggtitle("MCA plot with FactorMineR") + 
    scale_colour_discrete(name = "Variable") +
    geom_segment(data = mca_animal_exposed_end$vars_df, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
    theme_bw(20)+ xlab ("Dim 1 21.28 %")+ ylab("Dim 2 17.54%")
  

```

```{r animal-all-mca}
fviz_mca_var(mca_animal_all$mcares, col.var = 'cos2', alpha.var = 'contrib',
                           gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
                           select.var = list(cos2 = .6),
                           repel = TRUE, # Avoid text overlapping
                           ggtheme = theme_minimal()) 


arrows2 <- mca_animal_all$vars_df %>%
  filter(Variable %in% c('motility_adult', 'epibiotic', 'scavenger', 'filter_suspension', 'solitary_colonial')) %>%
  `row.names<-`(c('Colonial', 'Gregarious', 'Solitary', 'Motile','Sessile','Not Filter-Suspension', 'Filter or Suspension feeder',
                  'Not a Scavenger', 'Scavenger','Not Epibiotic','Epibiotic')) %>%
  filter(!(Variable == 'scavenger' & Dim.1 < 0 | Variable == 'epibiotic' & Dim.1 > 0 | Variable == 'filter_suspension' & Dim.1 > 0))



ggplot(data.frame(mca_animal_all$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = arrows2, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(arrows2), 
                                     colour = Variable), seed = 3824) +
  labs(title = "Animals: Sheltered 2010-present", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#94b4d8','#5484c4', '#374658', '#a08b5a', 'blue')) +
    geom_segment(data = arrows2, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.5,1.6) +
  ylim(-1.25,1) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220323_animal_all.png'), device = 'png', 
       height = 5, width = 7, units = 'in')  

```





```{r animal-mca-viz}
library(ggrepel)

#added scavenger and filter/suspension to >.6 cos2
# make subset of vars to plot arrows for
exposed_arrows1 <- mca_animal_exposed_start$vars_df %>%
  filter(Variable %in% c('motility_adult', 'epibiotic', 'filter_suspension')) %>%
  `row.names<-`(c('Motile','Sessile','Not Filter-Suspension', 'Filter or Suspension Feeder', 'Not Epibiotic','Epibiotic'))

exposed_arrows2 <- mca_animal_exposed_middle$vars_df %>%
  filter(Variable %in% c('motility_adult', 'epibiotic', 'filter_suspension', 'omnivore', 'predator', 'calcareous')) %>%
  `row.names<-`(c('Soft-bodied', 'Calcareous', 'Motile','Sessile', 'Not a predator', 'Predator', 'Not an omnivore', 'Omnivorous', 'Not Filter-Suspension','Filter or Suspension Feeder', 'Not Epibiotic','Epibiotic')) %>%
  filter(!(Variable == 'predator' & Dim.1 < 0 | Variable == 'calcareous' & Dim.1 < 0 | Variable == 'filter_suspension' & Dim.1 > 0))

exposed_arrows3 <- mca_animal_exposed_end$vars_df %>%
  filter(Variable %in% c('motility_adult', 'epibiotic', 'filter_suspension')) %>%
  `row.names<-`(c('Motile','Sessile','Not Filter-Suspension', 'Filter or Suspension Feeder', 'Not Epibiotic','Epibiotic'))

# exposed
ggplot(data.frame(mca_animal_exposed_start$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = exposed_arrows1, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(exposed_arrows1), 
                                     colour = Variable,
                                     family = 'quicksand'), seed = 3824) +
  labs(title = "Animals: 1982-1994, Exposed", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#94b4d8', '#5484c4')) +
    geom_segment(data = exposed_arrows1, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.5,1.6) +
  ylim(-1.25,1) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220327_animal_exposed_start.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
ggsave(filename = here(plots, '20220327_animal_exposed_start.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')

ggplot(data.frame(mca_animal_exposed_middle$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = exposed_arrows2, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(exposed_arrows2), 
                                     colour = Variable,
                                     family = 'quicksand'), seed = 3824) +
  labs(title = "Animals: 1996-2006, Exposed", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#543414','#806238', '#94b4d8', '#5484c4', '#9c8794', '#a08b5a', 'black', 'grey')) +
    geom_segment(data = exposed_arrows2, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.5,1.6) +
  ylim(-1.25,1.8) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220327_animal_exposed_middle.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
ggsave(filename = here(plots, '20220327_animal_exposed_middle.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')


ggplot(data.frame(mca_animal_exposed_end$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = exposed_arrows3, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(exposed_arrows3), 
                                     colour = Variable,
                                     family = 'quicksand'), seed = 3824) +
  labs(title = "Animals: 2010-present, Exposed", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#94b4d8','#5484c4')) +
    geom_segment(data = exposed_arrows3, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.5,1.6) +
  ylim(-1.25,1) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220327_animal_exposed_end.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
ggsave(filename = here(plots, '20220327_animal_exposed_end.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')







# make subset of vars to plot arrows for
sheltered_arrows1 <- mca_animal_sheltered_start$vars_df %>%
  filter(Variable %in% c('motility_adult', 'epibiotic')) %>%
  `row.names<-`(c('Motile','Sessile','Not Epibiotic','Epibiotic')) %>%
  filter(!(Variable == 'scavenger' & Dim.1 < 0 | Variable == 'epibiotic' & Dim.1 > 0 | Variable == 'filter_suspension' & Dim.1 > 0))

sheltered_arrows2 <- mca_animal_sheltered_middle$vars_df %>%
  filter(Variable %in% c('motility_adult', 'epibiotic', 'filter_suspension', 'predator', 'solitary_colonial')) %>%
  `row.names<-`(c('Colonial', 'Gregarious', 'Solitary', 'Motile','Sessile', 'Not a predator', 'Predator', 
                  'Not Filter-Suspension', 'Filter or Suspension feeder','Not Epibiotic','Epibiotic')) %>%
  filter(!(Variable == 'scavenger' & Dim.1 < 0 | Variable == 'epibiotic' & Dim.1 > 0 | 
             Variable == 'filter_suspension' & Dim.1 < 0 | Variable == 'filter_suspension' & Dim.1 > 0))

sheltered_arrows3 <- mca_animal_sheltered_end$vars_df %>%
  filter(Variable %in% c('motility_adult', 'epibiotic', 'scavenger', 'filter_suspension')) %>%
  `row.names<-`(c('Motile','Sessile','Not Filter-Suspension', 'Filter or Suspension feeder',
                  'Not a Scavenger', 'Scavenger','Not Epibiotic','Epibiotic')) %>%
  filter(!(Variable == 'scavenger' & Dim.1 < 0 | Variable == 'epibiotic' & Dim.1 > 0 | Variable == 'filter_suspension' & Dim.1 > 0))



# sheltered
ggplot(data.frame(mca_animal_sheltered_start$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = sheltered_arrows1, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(sheltered_arrows1), 
                                     colour = Variable,
                                     family = 'quicksand'), seed = 3824) +
  labs(title = "Animals: Sheltered 1982-1994", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#5484c4', 'green', 'orange', 'blue')) +
    geom_segment(data = sheltered_arrows1, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.5,1.6) +
  ylim(-1.25,1) +
  #scale_y_reverse(limits = c(1.3,-2.5)) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220327_animal_sheltered-start.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
ggsave(filename = here(plots, '20220327_animal_sheltered-start.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')


ggplot(data.frame(mca_animal_sheltered_middle$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = sheltered_arrows2, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(sheltered_arrows2), 
                                     colour = Variable,
                                     family = 'quicksand'), seed = 3824) +
  labs(title = "Animals: Sheltered 2010-present", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#5484c4', '#a08b5a', '#374658', 'orange', 'blue')) +
    geom_segment(data = sheltered_arrows2, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.5,1.6) +
  ylim(-1.25,1) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220327_animal_sheltered-middle.png'), device = 'png', 
       height = 5, width = 7, units = 'in')  
ggsave(filename = here(plots, '20220327_animal_sheltered-middle.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')  
  




ggplot(data.frame(mca_animal_sheltered_end$mcares$ind$coord), aes(x = Dim.1, y = Dim.2)) + 
    geom_hline(yintercept = 0, colour = "gray90") + 
    geom_vline(xintercept = 0, colour = "gray90") + 
    geom_point(colour = "gray50", alpha = 0.7) + 
    geom_density2d(colour = "gray87") + 
    geom_text_repel(data = sheltered_arrows3, aes(x = Dim.1, y = Dim.2, 
                                     label = rownames(sheltered_arrows3), 
                                     colour = Variable,
                                     family = 'quicksand'), seed = 3824) +
  labs(title = "Animals: Sheltered 2010-present", x = NULL, y = NULL) +
    scale_colour_manual(values = c('#806238','#94b4d8','#5484c4', '#3a230e', 'orange', 'blue')) +
    geom_segment(data = sheltered_arrows3, aes(x = 0, y = 0, xend = Dim.1, yend = Dim.2), 
               arrow = arrow(length = unit(0.2, "cm")), colour = "black") +
  xlim(-1.5,1.6) +
  ylim(-1.25,1) +
    theme_classic() +
  theme(axis.text  = element_blank(),
        legend.position = 'none')

ggsave(filename = here(plots, '20220327_animal_sheltered-end.png'), device = 'png', 
       height = 5, width = 7, units = 'in')  
ggsave(filename = here(plots, '20220327_animal_sheltered-end.pdf'), device = 'pdf', 
       height = 5, width = 7, units = 'in')  
  



as.data.frame(mca_animal_sheltered_end$mcares$ind$coord)


```

```{r animal-change}
ggplot() +
  geom_hline(yintercept = 0, colour = "gray90") + 
  geom_vline(xintercept = 0, colour = "gray90") + 
  geom_density2d(data = as.data.frame(mca_animal_sheltered_start$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "#94b4d8") +
  #geom_density2d(data = as.data.frame(mca_animal_sheltered_middle$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "grey60") +
  geom_density2d(data = as.data.frame(mca_animal_sheltered_end$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "#374658") +
  ylim(-1.2,1.3) +
  xlim(-1.5,1.5) +
  labs(title = 'Animals: Sheltered', subtitle = 'Beginning (blue) to end (black)', x = NA, y = NA) +
  theme_classic() +
  theme(axis.text = element_blank())

ggsave(filename = here(plots, '20220327_animal_sheltered_change.png'), device = 'png', 
       height = 5, width = 7, units = 'in')


ggplot() +
  geom_hline(yintercept = 0, colour = "gray90") + 
  geom_vline(xintercept = 0, colour = "gray90") + 
  geom_density2d(data = as.data.frame(mca_animal_exposed_start$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "#94b4d8") +
  #geom_density2d(data = as.data.frame(mca_animal_sheltered_middle$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "grey60") +
  geom_density2d(data = as.data.frame(mca_animal_exposed_end$mcares$ind$coord), aes(x = `Dim 1`, y = `Dim 2`), colour = "#374658") +
  ylim(-1.2,1.3) +
  xlim(-1.5,1.5) +
  labs(title = 'Animals: Exposed', subtitle = 'Beginning (blue) to end (black)', x = NA, y = NA) +
  theme_classic() +
  theme(axis.text = element_blank())

ggsave(filename = here(plots, '20220327_animal_exposed_change.png'), device = 'png', 
       height = 5, width = 7, units = 'in')
```


### other

```{r venns}
library(eulerr)

# make a list of lists
venn_all <- list(p_start = animal_sheltered_start$species,
                 p_end = animal_sheltered_end$species
                 #exposed = animal_exposed$species
                 )

# run eulerr
euleranimal <- euler(venn_all)

# plot output
plot(euleranimal,
     quantities = list(type = 'counts', cex = .75, fontface = 2),
     adjust_labels = TRUE,
     main = 'Unique fish species')
```