---
title: "Generating phys PC1"
author: "Ariana S Huffmyer, E5 RoL Team"
date: "2024"
output: 
  html_document:
    code_folding: hide
    toc: yes
    toc_depth: 6
    toc_float: yes
editor_options: 
  chunk_output_type: console
--- 

# Set Up    

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, warning=FALSE, message=FALSE)
```

```{r}
## install packages if you dont already have them
if (!require("tidyverse")) install.packages("tidyverse")
if (!require("ggplot2")) install.packages("ggplot2")
if (!require("RColorBrewer")) install.packages("RColorBrewer")
if (!require("lme4")) install.packages("lme4")
if (!require("lmerTest")) install.packages("lmerTest")
if (!require("car")) install.packages("car")
if (!require("effects")) install.packages("effects")
if (!require("ggfortify")) install.packages("ggfortify")
if (!require("cowplot")) install.packages("cowplot")
if (!require("vegan")) install.packages("vegan")
if (!require("corrr")) install.packages("corrr")
if (!require("ggcorrplot")) install.packages("ggcorrplot")
if (!require("GGally")) install.packages("GGally")
if (!require("broom")) install.packages("broom")
if (!require("cowplot")) install.packages("cowplot")
if (!require("pairwiseAdonis")) BiocManager::install("pairwiseAdonis")
#devtools::install_github("pmartinezarbizu/pairwiseAdonis/pairwiseAdonis")

# load packages
library(tidyverse)
library(ggplot2)
library(RColorBrewer)
library(lme4)
library(lmerTest)
library(car)
library(effects)
library(ggfortify)
library(cowplot)
library(vegan)
library(corrr)
library(ggcorrplot)
library(GGally)
library(broom)
library(cowplot)
library(pairwiseAdonis)
```

PICK OUT SAMPLES FOR MOLECULAR WORK 
GENERATE PC FOR EACH SPECIES 
GENERATE DATAFRAME OF PC VALUES FOR EACH SPECIES 

# Read in data 

```{r}
phys<-read_csv("data/phys_master_timeseries.csv")%>%
  select(colony_id_corr, species, timepoint, site, Host_AFDW.mg.cm2, Sym_AFDW.mg.cm2, Am, Rd, Ratio_AFDW.mg.cm2, cells.mgAFDW)%>%
  rename(colony_id=colony_id_corr)

length(unique(as.factor(phys$colony_id))) #140 colonies

metadata<-read_csv("data/rna_metadata.csv")%>%
  rename(colony_id=`Colony ID`, timepoint=Timepoint, haplotype=`Species/Strain`)%>%
  select(colony_id, timepoint, haplotype)%>%
  mutate(timepoint=gsub("TP", "timepoint", timepoint))
  
length(unique(as.factor(metadata$colony_id))) #30 colonies

data<-left_join(metadata, phys)

length(unique(as.factor(data$colony_id))) #30 colonies
```

Data now contains phys data for molecular samples.  

# Acropora PC1 

Generate a PC plot for physiology selected metrics for Acropora.  

```{r}
acr_data<-data%>%
  filter(species=="Acropora")

scaled_acr<-prcomp(acr_data%>%select(where(is.numeric)), scale=TRUE, center=TRUE) 

acr_PCA<-ggplot2::autoplot(scaled_acr, data=acr_data, frame.colour="timepoint", loadings=FALSE,  colour="timepoint", shape="site", loadings.label.colour="black", loadings.colour="black", loadings.label=TRUE, frame=FALSE, loadings.label.size=5, loadings.label.vjust=-1, size=5, frame.type = 'norm') + 
  theme_classic()+
   theme(legend.text = element_text(size=12), 
         legend.position="right",
        plot.background = element_blank(),
        legend.title = element_text(size=12, face="bold"), 
        axis.text = element_text(size=12), 
        axis.title = element_text(size=12,  face="bold"));acr_PCA

ggsave(acr_PCA, filename="output/02-phys-PC1/acr_PCA.png", width=6, height=6)
```

Extract PC scores for each sample on PC1. 

```{r}
acr_scores<-scaled_acr$x

acr_scores<-acr_scores%>%
  as.data.frame()%>%
  mutate(colony_id=acr_data$colony_id)%>%
  mutate(timepoint=acr_data$timepoint)%>%
  select(colony_id, timepoint, PC1)%>%
  write_csv(., file="output/02-phys-PC1/acr_PC1_scores.csv")
```

# Pocillopora PC1

Generate a PC plot for physiology selected metrics for Pocillopora  

```{r}
poc_data<-data%>%
  filter(species=="Pocillopora")

poc_data<-poc_data[complete.cases(poc_data), ]

scaled_poc<-prcomp(poc_data%>%select(where(is.numeric)), scale=TRUE, center=TRUE) 

poc_PCA<-ggplot2::autoplot(scaled_poc, data=poc_data, frame.colour="timepoint", loadings=FALSE,  colour="timepoint", shape="site", loadings.label.colour="black", loadings.colour="black", loadings.label=TRUE, frame=FALSE, loadings.label.size=5, loadings.label.vjust=-1, size=5, frame.type = 'norm') + 
  theme_classic()+
   theme(legend.text = element_text(size=12), 
         legend.position="right",
        plot.background = element_blank(),
        legend.title = element_text(size=12, face="bold"), 
        axis.text = element_text(size=12), 
        axis.title = element_text(size=12,  face="bold"));poc_PCA

ggsave(poc_PCA, filename="output/02-phys-PC1/poc_PCA.png", width=6, height=6)
```

Extract PC scores for each sample on PC1. 

```{r}
poc_scores<-scaled_poc$x

poc_scores<-poc_scores%>%
  as.data.frame()%>%
  mutate(colony_id=poc_data$colony_id)%>%
  mutate(timepoint=poc_data$timepoint)%>%
  select(colony_id, timepoint, PC1)%>%
  write_csv(., file="output/02-phys-PC1/poc_PC1_scores.csv")
```

# Porites PC1

Generate a PC plot for physiology selected metrics for Porites  

```{r}
por_data<-data%>%
  filter(species=="Porites")

por_data<-por_data[complete.cases(por_data), ]

scaled_por<-prcomp(por_data%>%select(where(is.numeric)), scale=TRUE, center=TRUE) 

por_PCA<-ggplot2::autoplot(scaled_por, data=por_data, frame.colour="timepoint", loadings=FALSE,  colour="timepoint", shape="site", loadings.label.colour="black", loadings.colour="black", loadings.label=TRUE, frame=FALSE, loadings.label.size=5, loadings.label.vjust=-1, size=5, frame.type = 'norm') + 
  theme_classic()+
   theme(legend.text = element_text(size=12), 
         legend.position="right",
        plot.background = element_blank(),
        legend.title = element_text(size=12, face="bold"), 
        axis.text = element_text(size=12), 
        axis.title = element_text(size=12,  face="bold"));por_PCA

ggsave(por_PCA, filename="output/02-phys-PC1/por_PCA.png", width=6, height=6)
```

Extract PC scores for each sample on PC1. 

```{r}
por_scores<-scaled_por$x

por_scores<-por_scores%>%
  as.data.frame()%>%
  mutate(colony_id=por_data$colony_id)%>%
  mutate(timepoint=por_data$timepoint)%>%
  select(colony_id, timepoint, PC1)%>%
  write_csv(., file="output/02-phys-PC1/por_PC1_scores.csv")
```
