##Blackspot for DRYAD library(dplyr) library(ggplot2) library(nlme) library(lme4) library(lubridate) ##set working directory to location of data files rm(list = ls()) setwd("C:/Users/Suzanne/Dropbox/Cody - black spot/dryad_data") #Read in the data file for blackspot infection of fish ##Info on data table: #FL is Fork Length of the fish, in mm ##Wt is weight of the fish, in g ##bs_infection is wether or not the fish was infected with black spot data <- read.csv("Steelhead_infections_data.csv") data$Date <- as.Date(data$Date, "%m/%d/%Y") data$juliandate <- julian(data$Date, origin = as.Date("2014-06-01")) ##Turn date into a numerical format for later use ##summarize number of fish infected or not infected by each sampling event fish_summary <- data %>% group_by(Site, site_name, Sample_event, Date, juliandate) %>% summarize(num_infected = sum(bs_infection == "Infected"), num_not = sum(bs_infection == "Not")) ##Read in data table for temperature data temp <- read.csv("stream_temperature_data.csv") temp$Date <- as.Date(temp$Date, "%m/%d/%Y") ##################################### ###Summarizing data into 7-day ADM### ##################################### ##first summarize by the maximum daily temperature tempsummary <- temp %>% group_by(Site, site_name, Date) %>% summarize(max_temp = max(temp)) ##next summarize the max temperature for the 7-day window associated with each fish sampling event fish_summary$ADM <- 0 ##create a new column to save the summary data ##Functions to summarize given a certain start date: early.dates.list <- function(date){c(date,date+1,date+2,date+3,date+4,date+5,date+6)} ##Early events are summarized by including the sample day and then following 6 days middle.dates.list <- function(date){c(date-3, date-2, date-1, date, date+1, date+2,date+3)} ##Middle events are summarized by centering the 7-day window on the sample day late.dates.list <- c(mdy("08/15/2014","08/16/2014","08/17/2014","08/18/2014","08/19/2014","08/20/2014","08/21/2014"))##Late events all occured on Aug 21, and include the 6 previous days ##For each line, this code subsets the temperature date for the early dates for the site of the sampling event, and then takes the average of those subsetted temperatures ##site1 site1 <- subset(tempsummary, Site == 1) fish_summary$ADM[1] <- mean(unlist(site1[(site1$Date %in% early.dates.list(fish_summary[1,4])),][,4])) fish_summary$ADM[2] <- mean(unlist(site1[(site1$Date %in% middle.dates.list(fish_summary[2,4])),][,4])) fish_summary$ADM[3] <- mean(unlist(site1[(site1$Date %in% late.dates.list),][,4])) #site2, datse are site2 <- subset(tempsummary, Site == 2) fish_summary$ADM[4] <- mean(unlist(site2[(site2$Date %in% early.dates.list(fish_summary[4,4])),][,4])) fish_summary$ADM[5] <- mean(unlist(site2[(site2$Date %in% middle.dates.list(fish_summary[5,4])),][,4])) fish_summary$ADM[6] <- mean(unlist(site2[(site2$Date %in% late.dates.list),][,4])) #site3 site3 <- subset(tempsummary, Site == 3) fish_summary$ADM[7] <- mean(unlist(site3[(site3$Date %in% early.dates.list(fish_summary[7,4])),][,4])) fish_summary$ADM[8] <- mean(unlist(site3[(site3$Date %in% middle.dates.list(fish_summary[8,4])),][,4])) fish_summary$ADM[9] <- mean(unlist(site3[(site3$Date %in% late.dates.list),][,4])) ##site 4 dates are 6/18/2014, 7/7/2014, 8/21/2014 site4 <- subset(tempsummary, Site == 4) fish_summary$ADM[10] <- mean(unlist(site4[(site4$Date %in% early.dates.list(fish_summary[10,4])),][,4])) fish_summary$ADM[11] <- mean(unlist(site4[(site4$Date %in% middle.dates.list(fish_summary[11,4])),][,4])) fish_summary$ADM[12] <- mean(unlist(site4[(site4$Date %in% late.dates.list),][,4])) #site5 site5 <- subset(tempsummary, Site == 5) fish_summary$ADM[13] <- mean(unlist(site5[(site5$Date %in% early.dates.list(fish_summary[13,4])),][,4])) fish_summary$ADM[14] <- mean(unlist(site5[(site5$Date %in% middle.dates.list(fish_summary[14,4])),][,4])) fish_summary$ADM[15] <- mean(unlist(site5[(site5$Date %in% late.dates.list),][,4])) #site6 site6 <- subset(tempsummary, Site == 6) fish_summary$ADM[16] <- mean(unlist(site6[(site6$Date %in% early.dates.list(fish_summary[16,4])),][,4])) fish_summary$ADM[17] <- mean(unlist(site6[(site6$Date %in% middle.dates.list(fish_summary[17,4])),][,4])) fish_summary$ADM[18] <- mean(unlist(site6[(site6$Date %in% late.dates.list),][,4])) #site7 site7 <- subset(tempsummary, Site == 7) fish_summary$ADM[19] <- mean(unlist(site7[(site7$Date %in% early.dates.list(fish_summary[19,4])),][,4])) fish_summary$ADM[20] <- mean(unlist(site7[(site7$Date %in% middle.dates.list(fish_summary[20,4])),][,4])) fish_summary$ADM[21] <- mean(unlist(site7[(site7$Date %in% late.dates.list),][,4])) rm(site1, site2, site3, site4,site5,site6,site7) ################################################################### ###___TEMPERATURE AND BLACK SPOT INFECTION RATES ON O. MYKISS___### ################################################################### ##Relating the 7-day ADM to infection rates ##Removing temperature observation that is artificially high because water level dropped, and logger was no longer submerged fish_summary$ADM[12]<- NA ##Mixed effects modelling between temperature, sample date (julian date), with Site as a random effect ##Model 1: Temperature as a fixed effect fit1 <- glmer(cbind(num_infected, num_not) ~ ADM + (1|Site), data = fish_summary, family = "binomial") summary(fit1) AIC(fit1) anova(fit1) ##Model 2: Sample date as fixed effect fit2 <- glmer(cbind(num_infected, num_not) ~ juliandate + (1|Site), data = fish_summary, family = "binomial") summary(fit2) AIC(fit2) ##Model 3: Temperature and sample date as fixed effects fit3 <- glmer(cbind(num_infected, num_not) ~ ADM + juliandate + (1|Site), data = fish_summary, family = "binomial") summary(fit3) AIC(fit3) ################################################################# ###____SNAIL COUNTS, TEMPERATURE, and O. MYKISS INFECTIONS___#### ################################################################# ##info on data table: ##AMM is Average Monthly Maximum, referring to the 30-day ADM, calculated from raw temperature data (available for download) snails <- read.csv("snail_count_data.csv") snails$Site <- as.factor(as.character(snails$Site)) snails$logsnails <- log(snails$num_snails+1) ##Exponential regression on snail abundance and 30-day ADM expfit <- lm(num_snails ~ exp(AMM), data = snails) summary(expfit) ##snails and O. mykiss infections ##summarizing full data set site_summary <- data %>% group_by(site_name) %>% summarize(num_infected = sum(bs_infection == "Infected"), num_not = sum(bs_infection == "Not")) snails <- merge(site_summary, snails, by = "site_name") ##relationship between snail abundance and steelhead infection rates (logistic regression) binomialfit <- glm(cbind(num_infected, num_not) ~ num_snails, data = snails, family = "binomial") summary(binomialfit) ############################################# ###_____BODY SIZE AND INFECTION RATES____#### ############################################# ##Subset the data for just the August sampling event, when we captured age-1+ fish in addition to age-0 fish augdat <- subset(data, Date == as.Date("8/21/2014", "%m/%d/%Y")) #two sided t-test between FL (fork length, or body size) and black spot infection status ttst <- t.test(data = augdat, FL ~ bs_infection, var.equal = F) ttst