Stats for Crab Responses to Predator Vocalizations

###May 13, 2014

###AR Hughes

### Stats for crab responses to predator vocalizations

install.packages(bbmle, dependencies=TRUE)

install.packages(lme4, dependencies=TRUE)

install.packages(lattice, dependencies=TRUE)

install.packages(multcomp, dependencies=TRUE)

install.packages(glmmADMB, dependencies=TRUE)

library(bbmle)

library(lme4)

library(lattice)

library(multcomp)

library(glmmADMB)

##########################

###Predator acoustic cue experiment###

##########################

###All days of experiment: Exposure day as random factor nested in Tank as random factor; week as random factor; Binomial distribution

acoustic = read.csv("acoustic_cue.csv", header=T) #read in file

names(acoustic)

acoustic$week=as.factor(acoustic$week) #convert to factor

acoustic$Tank=as.factor(acoustic$Tank)

acoustic$expos.day=as.factor(acoustic$expos.day)

model1 = glmmadmb(cbind(no..consumed, Total.clams.left)~1+(1|week)+(1|Tank/expos.day),family='binomial', data=acoustic)

model2 = glmmadmb(cbind(no..consumed, Total.clams.left)~Treatment+(1|week)+(1|Tank/expos.day),family='binomial', data=acoustic)

AICctab(logLik(model1),logLik(model2),delta=TRUE, weights=TRUE)

tukey=glht(model2, linfct=mcp(Treatment="Tukey")) #Tukey tests

summary.tukey = summary(tukey)

summary(tukey, test=adjusted("Westfall")) #Westfall adjustment of p-values for multiple comparisons

###Day 1 of experiment: Tank as random factor; week as random factor; Binomial distribution

acoustic_day1 = read.csv("acoustic_cue_day1only.csv", header=T) #read in file

names(acoustic_day1)

model1 = glmer(cbind(no..consumed, Total.clams.left)~1+(1|week)+(1|Tank), family=binomial, data=acoustic_day1)

model2 = glmer(cbind(no..consumed, Total.clams.left)~Treatment+(1|week)+(1|Tank), family=binomial, data=acoustic_day1)

AICctab(logLik(model1),logLik(model2),delta=TRUE, weights=TRUE)

##########################

###Predator acoustic cue vs. chemical cue experiment###

##########################

###Trial as random factor; Tank as random factor; Binomial distribution

acoustic_chemical = read.csv("acoustic_chemical.csv", header=T) #read in file

names(acoustic_chemical)

acoustic_chemical$trial=as.factor(acoustic_chemical$trial) #convert to factor

acoustic_chemical$tank=as.factor(acoustic_chemical$tank)

model1 = glmmadmb(cbind(Clams.eaten, Total.clams.left)~1+(1|trial)+(1|tank), family='binomial', data=acoustic_chemical)

model2 = glmmadmb(cbind(Clams.eaten, Total.clams.left)~sound+(1|trial)+(1|tank), family='binomial', data=acoustic_chemical)

model3 = glmmadmb(cbind(Clams.eaten, Total.clams.left)~water+(1|trial)+(1|tank), family='binomial', data=acoustic_chemical)

model4 = glmmadmb(cbind(Clams.eaten, Total.clams.left)~sound+water+(1|trial)+(1|tank), family='binomial', data=acoustic_chemical)

model5 = glmmadmb(cbind(Clams.eaten, Total.clams.left)~sound*water+(1|trial)+(1|tank),family='binomial', data=acoustic_chemical)

AICctab(logLik(model1),logLik(model2),logLik(model3),logLik(model4),logLik(model5),delta=TRUE, weights=TRUE)

acoustic_chemical$Inter=interaction(acoustic_chemical$sound, acoustic_chemical$water) #create an interaction variable

tukey=glht(model5, linfct=mcp(Inter="Tukey")) #Tukey tests

summary.tukey = summary(tukey)

summary(tukey, test=adjusted("Westfall")) # Westfall adjustment of p-values for multiple comparisons