-
Notifications
You must be signed in to change notification settings - Fork 0
/
CupRatio.R
213 lines (167 loc) · 9.65 KB
/
CupRatio.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
#Cup ratio & shell shape analysis and graphs
#RK 10/10/22
library(nlme)
library(lme4)
library(car)
library(MASS)
library(dplyr)
library(tidyverse)
library(ARTool)
library(plotrix)
library(lubridate)
library(hrbrthemes)
options(hrbrthemes.loadfonts = TRUE)
hrbrthemes::import_roboto_condensed()
library(patchwork)
#Import data
#allData<-read.csv("oysterDataAll.csv", na.strings=c(""," ","NA"))
allData<-read.csv("replicatetest.csv", na.strings=c(""," ","NA"))
#Fix date format
allData$Date<-mdy(allData$Date)
allData$Fan.ratio<-round(allData$Length/allData$Height, digits=2)
#Convert variables to factors
allData<-within(allData, {
Cage<-as.factor(Cage)
Bag<-as.factor(Bag)
Location<-as.factor(Location)
Gear<-as.factor(Gear)
Treatment<-as.factor(Treatment)
Replicate<-as.factor(Replicate)
})
str(allData)
allData$Replicate2<-paste0(allData$Treatment,".",allData$Replicate)
allData$Replicate2<-as.factor(allData$Replicate2)
str(allData)
SamplingFour<-allData[allData$Date=="2022-08-15",]
table(SamplingFour$Location,SamplingFour$Gear)
#DAY 1 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#Select data
SamplingOne<-allData[allData$Date=="2022-06-14",]
str(SamplingOne)
table(SamplingOne$Location, SamplingOne$Gear)
#Combined graph - cup ratio
cupRatioGraph1<-ggplot(data = SamplingOne, aes(x = Gear, y = Cup.ratio, fill=Location))+geom_boxplot()+ylab("Cup ratio (shell width/height)")+theme_classic()#+scale_y_continuous(limits=c(0,0.5))
#Combined graph - shell shape
shellShapeGraph1<-ggplot(data = SamplingOne, aes(x = Gear, y = Shell.shape, fill=Location))+geom_boxplot()+ylab("Shell shape")+theme_classic()+scale_y_continuous(limits=c(0,9))
#Graph inside vs. graph outside - cup ratio
cupRatioGraph2 <- ggplot(SamplingOne, aes(x=Gear, y=Cup.ratio, group=Gear)) +
geom_boxplot(aes(fill=Gear))
cupRatioGraph2 + facet_grid(. ~ Location)
#Graph inside vs. graph outside - shell shape
shellShapeGraph2 <- ggplot(SamplingOne, aes(x=Gear, y=Shell.shape, group=Gear)) +
geom_boxplot(aes(fill=Gear))
shellShapeGraph2 + facet_grid(. ~ Location)
#Cup ratio ART - significant effect of both location and gear
artCupRatioOne<-art(Cup.ratio ~ Gear * Location, data=SamplingOne)
artCupRatioOne #appropriate
anova(artCupRatioOne)
#Shell shape ART - significant effect of location
artShellShapeOne<-art(Shell.shape ~ Gear * Location, data=SamplingOne)
artShellShapeOne #appropriate
anova(artShellShapeOne)
#Location cup ratio post-hoc - inside higher cup ratio than outside
LocationPostHoc<-art.con(artCupRatioOne, "Location", adjust="bonferroni")# %>% run ART-C for X1 × X2
summary(LocationPostHoc) %>% #add significance stars to the output
mutate(sig. = symnum(p.value, corr=FALSE, na=FALSE,
cutpoints = c(0, .001, .01, .05, .10, 1),
symbols = c("***", "**", "*", ".", " ")))
#Gear post-hoc - FC cup ratio higher than BP
GearPostHoc<-art.con(artCupRatioOne, "Gear", adjust="bonferroni")# %>% run ART-C for X1 × X2
summary(GearPostHoc) %>% #add significance stars to the output
mutate(sig. = symnum(p.value, corr=FALSE, na=FALSE,
cutpoints = c(0, .001, .01, .05, .10, 1),
symbols = c("***", "**", "*", ".", " ")))
#Location shell shape post-hoc - inside
LocationPostHocShell<-art.con(artShellShapeOne, "Location", adjust="bonferroni")# %>% run ART-C for X1 × X2
summary(LocationPostHocShell) %>% #add significance stars to the output
mutate(sig. = symnum(p.value, corr=FALSE, na=FALSE,
cutpoints = c(0, .001, .01, .05, .10, 1),
symbols = c("***", "**", "*", ".", " ")))
#DAY 2 Cup Ratio ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#Import and select data
SamplingTwo<-allData[allData$Date=="2022-07-05",]
SamplingThree<-allData[allData$Date=="2022-07-26",]
#Summary function
data_summary <- function(data, varname, groupnames){
summary_func <- function(x, col){
c(mean = mean(x[[col]], na.rm=TRUE),
SE = std.error(x[[col]], na.rm=TRUE))
}
data_sum<-plyr::ddply(data, groupnames, .fun=summary_func,
varname)
return(data_sum)
}
#Cup ratio time graph
timeGraphCupRatioDf<-data_summary(allData, "Cup.ratio",
groupnames=c("Date", "Location", "Gear"))
timeGraphCupRatio<-ggplot(timeGraphCupRatioDf, aes(x=Date, y=mean, color=Location, linetype=Gear)) + geom_line() +geom_point()+theme_classic()+geom_errorbar(aes(ymin=mean-SE, ymax=mean+SE), width=.2, position=position_dodge(0.05))
timeGraphCupRatio
#Shell shape time graph
timeGraphShellShapeDf<-data_summary(allData, "Shell.shape",groupnames=c("Date", "Location", "Gear"))
timeGraphShellShape<-ggplot(timeGraphShellShapeDf, aes(x=Date, y=mean, color=Location, linetype=Gear)) +geom_line() +geom_point()+theme_classic()+geom_errorbar(aes(ymin=mean-SE, ymax=mean+SE), width=.2,position=position_dodge(0.05))+ylab("Shell shape index")+xlab("")
timeGraphShellShape
#Fan ratio time graph
timeGraphFanRatioDf<-data_summary(allData, "Fan.ratio",groupnames=c("Date", "Location", "Gear"))
timeGraphFanRatio<-ggplot(timeGraphFanRatioDf, aes(x=Date, y=mean, color=Location, linetype=Gear)) + geom_line() +geom_point()+theme_classic()+geom_errorbar(aes(ymin=mean-SE, ymax=mean+SE), width=.2,position=position_dodge(0.05))
timeGraphFanRatio
#Day one fan ratio - both location and gear significant
artFanRatio<-art(Cup.ratio ~ Gear * Location + (1|Replicate2), data=SamplingOne)
artFanRatio #not appropriate
anova(artFanRatio)
fanRatioGraph1<-ggplot(data = SamplingOne, aes(x = Gear, y = Fan.ratio, fill=Location))+geom_boxplot()+scale_y_continuous(limits=c(0,1.5))+ylab("Fan ratio (shell length/height)")+theme_classic()
#DAY 4 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
SamplingFour<-allData[allData$Date=="2022-08-15",]
#Combined graph cup ratio
ClarkFEST1<-ggplot(data = SamplingFour, aes(x = Gear, y = Cup.ratio, fill=Location))+geom_boxplot()+ylab("Cup ratio (SW/SH)")+theme_classic()+theme(axis.title.y = element_text(margin = margin(r = 15)))
ClarkFEST2<-ggplot(data = SamplingFour, aes(x = Gear, y = Shell.shape, fill=Location))+geom_boxplot()+ylab("Shell shape index")+theme_classic()+ theme(axis.title.y = element_text(margin = margin(r = 15)))
ClarkFEST2
ClarkFEST<- ClarkFEST1 + ClarkFEST2+ plot_layout(nrow=1, guides = "collect")+ plot_annotation(tag_levels = 'A') & theme(legend.position = "bottom")
ClarkFEST
ggplot(data = SamplingFour, aes(x = Replicate2, y = Cup.ratio, fill=Treatment))+geom_boxplot()+ylab("Cup ratio (shell width/height)")+theme_classic()#+scale_y_continuous(limits=c(0,0.5))
#Cup ratio ART 4- significant effect of gear, location, and interaction
artCupRatioFour<-art(Cup.ratio ~ Gear * Location + (1|Replicate2), data=SamplingFour)
artCupRatioFour #appropriate
anova(artCupRatioFour)
test2<-data_summary(SamplingFour, "Cup.ratio", groupnames=c("Treatment"))
test2
test1<-data_summary(allData, "Cup.ratio", groupnames=c("Treatment"))
test1
#Combined graph cup ratio
ggplot(data = allData, aes(x = Gear, y = Cup.ratio, fill=Location))+geom_boxplot()+ylab("Cup ratio (shell width/height)")+theme_classic()#+scale_y_continuous(limits=c(0,0.5))
artCupRatioAll<-art(Cup.ratio ~ Gear * Location + (1|Date), data=allData)
artCupRatioAll #appropriate
anova(artCupRatioAll)
#Gear post-hoc - BP > FB
CupAllGearPostHoc<-art.con(artCupRatioAll, "Gear", adjust="bonferroni")
summary(CupAllGearPostHoc) %>% #add significance stars to the output
mutate(sig. = symnum(p.value, corr=FALSE, na=FALSE,
cutpoints = c(0, .001, .01, .05, .10, 1),
symbols = c("***", "**", "*", ".", " ")))
#interaction post-hoc: BPi>FBi, BPi>FCi, FCo>FBi
CupAllGearPostHoc<-art.con(artCupRatioAll, "Gear:Location", adjust="bonferroni")
summary(CupAllGearPostHoc) %>% #add significance stars to the output
mutate(sig. = symnum(p.value, corr=FALSE, na=FALSE,
cutpoints = c(0, .001, .01, .05, .10, 1),
symbols = c("***", "**", "*", ".", " ")))
#Shell shape ART four- significant effect of gear
artShellFour<-art(Shell.shape ~ Gear * Location, data=SamplingFour)
anova(artShellFour)
#floating gear higher than BP
Shell3PostHoc<-art.con(artShellThree, "Gear", adjust="bonferroni")
summary(Shell3PostHoc) %>% #add significance stars to the output
mutate(sig. = symnum(p.value, corr=FALSE, na=FALSE,
cutpoints = c(0, .001, .01, .05, .10, 1),
symbols = c("***", "**", "*", ".", " ")))
# initialFBi.1<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FBi.1","Cup.ratio"])
# initialFBi.2<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FBi.2","Cup.ratio"])
# initialFBi.3<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FBi.3","Cup.ratio"])
# initialFBo.1<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FBo.1","Cup.ratio"])
# initialFBo.2<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FBo.2","Cup.ratio"])
# initialFBo.3<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FBo.3","Cup.ratio"])
#
# initialFCi.1<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FCi.1","Cup.ratio"])
# initialFCi.2<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FCi.2","Cup.ratio"])
# initialFCi.3<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FCi.3","Cup.ratio"])
# initialFCo.1<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FCo.1","Cup.ratio"])
# initialFCo.2<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FCo.2","Cup.ratio"])
# initialFCo.3<-mean(allData[allData$Date=="2022-06-14" &allData$Replicate2=="FCo.3","Cup.ratio"])