exec/RepeatedMeasuresAnova_forRoots.R

In bmcnellis/SDEF.analysis: Data Analysis for 2016 San Dimas Experimental Forest Project

#### Repeated measures anova on root biomass data
library(ggplot2)
library(gridExtra)
library(cowplot)
library(lme4)
library(lmerTest)
library(lubridate)
library(emmeans)
setwd("D:/R_Workspace/SDEF/SDEF.analysis/data")
#load('SDEF_roots_post_aggregation_by_month.rda')
load("Alldates_depthbin_rootbiomass.Rdata")
load("SDEF_roots_pre_aggregation.rda")
###
biodensity <- SDEF_roots_pre_aggregation
#biodensity$mean_biodensity_dm <- biodensity$mean_biodensity*1000
#biodensity$sd_biodensity_dm <- biodensity$sd_biodensity*1000
biodensity$date <- as.POSIXct(biodensity$date, format = "%m/%d/%Y")
biodensity$month <- as.factor(month(biodensity$date))
biodensity$depth_bin <- as.factor(biodensity$depth_bin)
biodensity[1826, ] <- biodensity[ 2, ]
biodensity[1827, ] <- biodensity[ 2, ]
biodensity[1828, ] <- biodensity[ 2, ]
biodensity[c(1826:1828), 6] <- 1
biodensity[c(1826:1828), c(4,5,7)] <- 0
biodensity[1827, 2] <- 2
biodensity[1828, 2] <- 3
biodensity <- subset(biodensity, biodensity$month != 10)

# repeated measures  ANOVA Biomass
biomod1 <- lmer(length ~ site + depth_bin + month +  site*month*depth_bin + (1|tube) +(1|date), data = biodensity)
step(biomod1)
plot(biomod1)
anova(biomod1)

emmip(biomod1, site ~ depth_bin)
emout <- emmeans(biomod1, list(pairwise ~ site|depth_bin|month), adjust = "Tukey")
emout1 <- as.data.frame(emout$`pairwise differences of site | depth_bin, month`)

#write.csv(emout1, file = "eemeans_roots_veg_month_depthoutput_04_2019.csv")
###just depth
depth_mod <- lmer(diameter ~ depth_bin + site + site*depth_bin + (1|tube) +(1|date), data = biodensity)
plot(depth_mod)
anova(depth_mod)
emmeans(depth_mod, list(pairwise ~ site|depth_bin), adjust = "Tukey")
### Root length  and diameter
biomod1 <- lmer(length ~ site + depth_bin + month +  site*month + site:depth_bin+ (1|tube) +(1|date), data = biodensity)
plot(biomod1)
anova(biomod1)

emmip(biomod1, site ~ depth_bin)
emmeans(biomod1, list(pairwise ~ site|month), adjust = "Tukey")
emout1 <- as.data.frame(emout$`pairwise differences of site | depth_bin, month`)

#write.csv(emout1, file = "eemeans_roots_length_veg_month_depthoutput_04_2019.csv")

# try edies suggestion with separate anovas for each depth
dep1 <- subset(biodensity, biodensity$depth_bin == 1)
mod <- lmer(length ~ site +  month +  site*month + (1|tube) +(1|date), data = dep1)
anova(mod)
plot(mod)
emmeans(mod, list(pairwise ~ site|month),  adjust = "Tukey")

dep2 <- subset(biodensity, biodensity$depth_bin == 2)
mod <- lmer(length ~ site +  month +  site*month + (1|tube) +(1|date), data = dep2)
anova(mod)
plot(mod)
emmeans(mod, list(pairwise ~ site|month),  adjust = "Tukey")

dep3 <- subset(biodensity, biodensity$depth_bin == 3)
mod <- lmer(length ~ site +  month +  site*month + (1|tube) +(1|date), data = dep3)
anova(mod)
plot(mod)

dep4 <- subset(biodensity, biodensity$depth_bin == 4)
mod <- lmer(length ~ site +  month +  site*month + (1|tube) +(1|date), data = dep4)
anova(mod)
plot(mod)
emmeans(mod, list(pairwise ~ site|month),  adjust = "Tukey")


#### Review #2 Suggestion
biodensity$jday <- lubridate::yday(biodensity$date)

mod1 <- lmer(length ~ site + jday + site*jday +(1|tube)+ (1|date), data = dep1)
anova(mod1)
plot(mod1)