Set working directory: Please paste the whole path to the folder where you have your files
setwd("C:/Users/MariaNelly/Documents/Dr_Mercer_Project/01_CRD")
1 Read your data set and double check that dependent and indepent variables are correctly read by R
CRD_Data <- data.frame(read.csv("CRD_Data.csv",header=T))
str(CRD_Data)
## 'data.frame': 9 obs. of 4 variables:
## $ Plot : int 1 2 3 4 5 6 7 8 9
## $ Treatment : Factor w/ 3 levels "Ambient moisture",..: 3 3 1 2 1 1 3 2 2
## $ Replication: int 1 3 1 2 3 2 2 1 3
## $ biomass : num 121 118 111 131 105 ...
Perform ANOVA for CRD
Anova_CRD <- aov(biomass~Treatment,data=CRD_Data)
summary(Anova_CRD)
Df Sum Sq Mean Sq F value Pr(>F)
Treatment 2 1058.0 529.0 58.04 0.000119 ***
Residuals 6 54.7 9.1
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Run a Tukey test
TukeyHSD(Anova_CRD)
## Tukey multiple comparisons of means
## 95% family-wise confidence level
##
## Fit: aov(formula = biomass ~ Treatment, data = CRD_Data)
##
## $Treatment
## diff lwr upr p adj
## Heavy Irrigation-Ambient moisture 26.50000 18.936659 34.063341 0.0000944
## Light irrigation-Ambient moisture 14.76667 7.203326 22.330008 0.0023563
## Light irrigation-Heavy Irrigation -11.73333 -19.296674 -4.169992 0.0074656
Run an LSD test
library(agricolae)
## Warning: package 'agricolae' was built under R version 3.3.3
LSD_out <- LSD.test(Anova_CRD,"Treatment",alpha=0.05,p.adj="bonferroni")
LSD_out
## $statistics
## Mean CV MSerror LSD
## 120.2556 2.510498 9.114444 8.103625
##
## $parameters
## Df ntr bonferroni alpha test name.t
## 6 3 3.287455 0.05 bonferroni Treatment
##
## $means
## biomass std r LCL UCL Min Max
## Ambient moisture 106.5000 3.675595 3 102.2350 110.7650 103.5 110.6
## Heavy Irrigation 133.0000 1.915724 3 128.7350 137.2650 130.8 134.3
## Light irrigation 121.2667 3.187998 3 117.0016 125.5317 118.4 124.7
##
## $comparison
## NULL
##
## $groups
## trt means M
## 1 Heavy Irrigation 133.0000 a
## 2 Light irrigation 121.2667 b
## 3 Ambient moisture 106.5000 c
Draw a preliminary boxplot by treatment
boxplot(biomass~Treatment, data=CRD_Data)
Use ggplot2 to make a boxplot
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 3.3.3
CRD_plot <- ggplot(CRD_Data,aes(x=Treatment,y=biomass))+geom_boxplot()
CRD_plot + scale_y_continuous (name="Biomass in Kilograms") + theme_bw()