6_modeling_analysis.Rmd

---
title: "Modeling analysis of E5 time series biological data"
author: "Ariana S Huffmyer, E5 RoL Team"
date: "03/21/2023"
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("Hmisc")) install.packages("Hmisc")

# 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(Hmisc)
```

# Load dataframe

Load in master dataframe generated from 1_assemble_data.Rmd.  
```{r}
master<-read.csv("Output/master_timeseries.csv")
```

# Run models 

## All species 

Model the influence of physiological metrics on calcification. 
```{r}
df<-master%>%
  select(Host_AFDW.mg.cm2, Sym_AFDW.mg.cm2, Am, AQY, Rd, calc.umol.cm2.hr, cells.mgAFDW, prot_mg.mgafdw, cre.umol.mgafdw, Ratio_AFDW.mg.cm2, Total_Chl_cell, Total_Chl, Ik, Ic)

model1<-lm(calc.umol.cm2.hr~ ., data=df)

summary(model1)
```


Show effect plot
```{r}
plot(allEffects(model1)) 
```

The effect plot show that calcification decreases with increasing host biomass, increases with increasing max photosynthesis, and increases with increasing antioxidant capacity. The surprising negative relationship between host biomass and calcification could be driven by Porites, which is a species with high biomass. We will separate out species below. 

Next, look for covariates and simplify the model by removing highly correlated variables.  

```{r}
cor_all<-rcorr(as.matrix(df))
#combine all values 
flattenCorrMatrix <- function(cormat, pmat) {
  ut <- upper.tri(cormat)
  data.frame(
    row = rownames(cormat)[row(cormat)[ut]],
    column = rownames(cormat)[col(cormat)[ut]],
    cor  =(cormat)[ut],
    p = pmat[ut]
    )
}

cor_all_df<-flattenCorrMatrix(cor_all$r, cor_all$P)
```

List correlations with absolute value r>0.8.  
```{r}
cor_all_df%>%
  filter(abs(cor)>0.80)
```
Respiration and Am are > 0.8.  

Re run the model without respiration (not significant in original model). 
```{r}
df<-master%>%
  select(Host_AFDW.mg.cm2, Sym_AFDW.mg.cm2, Am, AQY, Ik, Ic, calc.umol.cm2.hr, cells.mgAFDW, prot_mg.mgafdw, cre.umol.mgafdw, Ratio_AFDW.mg.cm2, Total_Chl_cell, Total_Chl)

model1<-lm(calc.umol.cm2.hr~ ., data=df)

summary(model1)
```

The results are the same as the above model. 

List correlations with high significance p<0.01.  
```{r}
cor_all_df%>%
  filter(p<0.01)
```
Many correlations are highly significant.  

View the correlations between calcification and other physiological parameters.  
```{r}
cor_all_df%>%
  filter(p<0.01)%>%
  filter(row=="calc.umol.cm2.hr")
```
There are significant, but weak positive correlations exist between calcification and host protein, antioxidant capacity, and total chlorophyll (negative). 

## Acropora 

Model the influence of physiological metrics on calcification. 
```{r}
df_acr<-master%>%
  filter(species=="Acropora")%>%
  select(Host_AFDW.mg.cm2, Sym_AFDW.mg.cm2, Am, AQY, Rd, Ik, Ic, calc.umol.cm2.hr, cells.mgAFDW, prot_mg.mgafdw, cre.umol.mgafdw, Ratio_AFDW.mg.cm2, Total_Chl_cell, Total_Chl)

model1_acr<-lm(calc.umol.cm2.hr~ ., data=df_acr)

summary(model1_acr)
```

Show effect plot
```{r}
plot(allEffects(model1_acr)) 
```

There is a positive relationship with calcification and cell density, host protein, Ik, Ic, AQY, and chlorophyll per cell. Interestingly, there is a negative relationship between total chlorophyll per mg afdw, Am, and calcification.  

Next, look for covariates and simplify the model by removing highly correlated variables.  

```{r}
cor_acr<-rcorr(as.matrix(df_acr))
#combine all values 

cor_acr_df<-flattenCorrMatrix(cor_acr$r, cor_acr$P)
```

List correlations with absolute value r>0.9.  
```{r}
cor_acr_df%>%
  filter(abs(cor)>0.8)
```
No correlations > 0.8. The model will not be revised.  

List correlations with high significance p<0.01.  
```{r}
cor_acr_df%>%
  filter(p<0.01)
```
Many correlations are highly significant. 

View the correlations between calcification and other physiological parameters.  
```{r}
cor_acr_df%>%
  filter(p<0.01)%>%
  filter(row=="calc.umol.cm2.hr")
```
There are significant positive correlations between calcification and protein and S:H biomass ratio. 

## Pocillopora 

Model the influence of physiological metrics on calcification. 
```{r}
df_poc<-master%>%
  filter(species=="Pocillopora")%>%
  select(Host_AFDW.mg.cm2, Sym_AFDW.mg.cm2, Am, AQY, Rd, Ik, Ic, calc.umol.cm2.hr, cells.mgAFDW, prot_mg.mgafdw, cre.umol.mgafdw, Ratio_AFDW.mg.cm2, Total_Chl_cell, Total_Chl)

model1_poc<-lm(calc.umol.cm2.hr~ ., data=df_poc)

summary(model1_poc)
```

Antioxidant capacity is the only significant predictor of calcification.  

Show effect plot
```{r}
plot(allEffects(model1_poc)) 
```

There is a negative relationship between antioxidant capacity and calcification in Pocillopora.  

Next, look for covariates and simplify the model by removing highly correlated variables.  

```{r}
cor_poc<-rcorr(as.matrix(df_poc))
#combine all values 

cor_poc_df<-flattenCorrMatrix(cor_poc$r, cor_poc$P)
```

List correlations with absolute value r>0.8.  
```{r}
cor_poc_df%>%
  filter(abs(cor)>0.8)
```
No correlations > 0.8. The model will not be adjusted.   

List correlations with high significance p<0.01.  
```{r}
cor_poc_df%>%
  filter(p<0.01)
```
Many correlations are highly significant. 

View the correlations between calcification and other physiological parameters.  
```{r}
cor_poc_df%>%
  filter(p<0.01)%>%
  filter(row=="calc.umol.cm2.hr")
```
There is one significant negative correlation between calcification and antioxidant capacity.   

## Porites 

Model the influence of physiological metrics on calcification. 
```{r}
df_por<-master%>%
  filter(species=="Porites")%>%
  select(Host_AFDW.mg.cm2, Sym_AFDW.mg.cm2, Am, AQY, Rd, Ik, Ic, calc.umol.cm2.hr, cells.mgAFDW, prot_mg.mgafdw, cre.umol.mgafdw, Ratio_AFDW.mg.cm2, Total_Chl_cell, Total_Chl)%>%
  filter()

model1_por<-lm(calc.umol.cm2.hr~ ., data=df_por)

summary(model1_por)
```


Show effect plot
```{r}
plot(allEffects(model1_poc)) 
```


Next, look for covariates and simplify the model by removing highly correlated variables.  

```{r}
cor_por<-rcorr(as.matrix(df_por))
#combine all values 

cor_por_df<-flattenCorrMatrix(cor_por$r, cor_por$P)
```

List correlations with absolute value r>0.8.  
```{r}
cor_por_df%>%
  filter(abs(cor)>0.8)
```
No correlations > 0.8. The model will not be revised.   

List correlations with high significance p<0.01.  
```{r}
cor_por_df%>%
  filter(p<0.01)
```
Many correlations are highly significant. 

View the correlations between calcification and other physiological parameters.  
```{r}
cor_por_df%>%
  filter(p<0.01)%>%
  filter(row=="calc.umol.cm2.hr")
```
There is nothing that correlates with calcification in Porites.