####DROSOPHILA GLUE PREVENTS FROM PREDATION #************************************************************************************************************************************************************************************* #************************************************************************************************************************************************************************************* ########################################################ADHESION OF VINCENNES SPECIES##################### rm(list=ls(all=T)) #install.packages("ggplot2") #install.packages("gridExtra") library(ggplot2) library('gridExtra') ############################ #total_table <- read.csv("./data/To_Dryad/pupa_adhesion_print.csv", sep=",", fileEncoding = "UTF-7") total_table <- read.table("pupa_adhesion_print.txt", sep='\t', header = TRUE) #Only select "ok" assays total_table_ok <- subset(total_table, total_table$Comment_on_this_sample == "ok") ################################################ #Plot according to species total_table_ok$Stock <- factor(total_table_ok$Stock, levels = c('suzukii_vincennes','simulans_vincennes', 'hydei_vincennes',ordered = TRUE)) n_fun <- function(total_table_ok){ return(data.frame(y = 1100, label = paste0("n = ",length(total_table_ok)))) } p1 <- ggplot2::ggplot(data= total_table_ok, aes(x= Stock, y= force_detachment_mN)) p1 <- p1 + geom_boxplot(width= 0.4, outlier.colour = "white", fill = "lightgrey") + xlab("Stock") + ylab("Force_mN") p1 <- p1 + geom_point(alpha=0.4,shape = 21, colour = "black", size = 2, stroke = 1, position=position_jitter(0.02)) #p1 <- p1 + stat_summary(fun.data = n_fun, geom = "text", size=5) p1 <- p1 + theme_classic() + theme(plot.title = element_text(color="black", size=20, face = "bold"), axis.title = element_text(size = 20), axis.text = element_text(color="black", size=15)) p1 <- p1 + scale_y_continuous(name= "Force (mN)", breaks=seq(0,1400,200)) + xlab("Species") p1 #Median values #library(asbio) #install.packages("asbio") median(subset(total_table_ok, total_table_ok$Stock == 'hydei_vincennes')$force_detachment_mN, na.rm =T) median(subset(total_table_ok, total_table_ok$Stock == 'simulans_vincennes')$force_detachment_mN, na.rm =T) mad(subset(total_table_ok, total_table_ok$Stock == 'simulans_vincennes')$force_detachment_mN, na.rm =T) median(subset(total_table_ok, total_table_ok$Stock == 'suzukii_vincennes')$force_detachment_mN, na.rm =T) mad(subset(total_table_ok, total_table_ok$Stock == 'suzukii_vincennes')$force_detachment_mN, na.rm =T) #test normality shapiro.test(total_table_ok$force_detachment_mN) shapiro.test(log(total_table_ok$force_detachment_mN[resultat_force$Stock == 'hydei_vincennes'])) shapiro.test(log(total_table_ok$force_detachment_mN[resultat_force$Stock == 'simulans_vincennes'])) shapiro.test(log(total_table_ok$force_detachment_mN[resultat_force$Stock == 'suzukii_vincennes'])) # Not all species are normally distributed suz <- subset(total_table_ok,total_table_ok$Stock == "suzukii_vincennes") #Distribution of forces p4 <- ggplot(data= total_table_ok, aes(force_detachment_mN)) + geom_histogram() p4 <- p4 + facet_wrap(~ Stock) p4 #As it is not normal: test if Force depends on species kruskal.test(force_detachment_mN ~ Stock, data= total_table_ok) pairwise.wilcox.test(total_table_ok$force_detachment_mN, total_table_ok$Stock, p.adjust.method ="holm") ######################################################################### #Correlation Force VS Pupa-substrate contact area total_table_ok$Area_px <- as.numeric(as.character(total_table_ok$Area_px)) total_table_ok$Area_mm <- (total_table_ok$Area_px)*(1/(total_table_ok$scale_px/total_table_ok$scale_mm)^2) #Plot Force VS pupa-substrate contact p2 <- ggplot2::ggplot(data= total_table_ok, aes(x= Area_mm, y= force_detachment_mN)) p2 <- p2 + geom_point(aes(shape=Stock), size=2) p2 <- p2 + theme_classic() + theme(plot.title = element_text(color="black", size=20, face = "bold"), axis.title = element_text(size = 20), axis.text = element_text(color="black", size=15)) p2 <- p2 + scale_shape_manual(values=c(1, 17, 8)) + scale_y_continuous(name= "Force (mN)", breaks=seq(0,1400,200)) p2 <- p2 + theme(legend.position = c(0.7, 0.1)) p2 ################################################################ ### Correct Force by contact area n_fun <- function(total_table_ok){ return(data.frame(y = 1100, label = paste0("n = ",length(total_table_ok)))) } p3 <- ggplot2::ggplot(data= total_table_ok, aes(x= Stock, y= force_detachment_mN/Area_mm)) p3 <- p3 + geom_boxplot(width= 0.4, outlier.colour = "white", fill = "lightgrey") + xlab("Stock") + ylab("Force_mN") p3 <- p3 + geom_point(alpha=0.4,shape = 21, colour = "black", size = 2, stroke = 1, position=position_jitter(0.02)) p3 <- p3 + stat_summary(fun.data = n_fun, geom = "text", size=5) p3 <- p3 + theme_classic() + theme(plot.title = element_text(color="black", size=20, face = "bold"), axis.title = element_text(size = 20), axis.text = element_text(color="black", size=15)) p3 <- p3 + scale_y_continuous(name= "Force (mN)", breaks=seq(0,2100,200)) + xlab("Species") p3 #test normality shapiro.test((total_table_ok$force_detachment_mN/total_table_ok$Area_mm)[total_table_ok$Stock == 'hydei_vincennes']) shapiro.test((total_table_ok$force_detachment_mN/total_table_ok$Area_mm)[total_table_ok$Stock == 'suzukii_vincennes']) shapiro.test((total_table_ok$force_detachment_mN/total_table_ok$Area_mm)[total_table_ok$Stock == 'simulans_vincennes']) # Not all species are normally distributed #As it is not normal: test if Force depends on species kruskal.test(force_detachment_mN/Area_mm ~ Stock, data= total_table_ok) pairwise.wilcox.test(total_table_ok$force_detachment_mN/total_table_ok$Area_mm, total_table_ok$Stock, p.adjust.method ="holm") ###########FIGURE 1 grid.arrange(p1,p2,p3, ncol=3) ########### ######################################################################################################## #REGRESSION FORCES VS AREAS #Use sma function to estimate correlations library(smatr) ### Several sample analyses ### # Fit SMA's separately at each of Stocks, # and test for common slope: com.test <- sma(force_detachment_mN~Area_mm*Stock, data=total_table_ok) com.test #slopes are equal #summary() gives Correlation R squared and p-value and elevation and slope for each group summary(com.test) #correlations are significative for each group # Plot Force VS Area separately for each group: plot(com.test) #Estimate SMA common slope estimate with(total_table_ok, slope.com(force_detachment_mN,Area_mm, Stock, method='SMA')) #b = 491 # Fit SMA's separately at each of the stocks and test for common elevation: sma(force_detachment_mN~Area_mm+Stock, data=total_table_ok) #elevation is different #Produce a residual plot and qq-plot to check assumptions (residuals with normal distribution): plot(com.test, which = 'res') #Produce a normal quantil plot plot(com.test, which="qq") #Residuals seem to follow normal distribution ######################################################################################################### #RATIO AREA PRINT / AREA PUPA pupal_shape <- read.table("pupa_adhesion_size.txt", sep="\t", header=T) #Convert area_px int square mm #Convert px to squared mm pupal_shape$area <- as.numeric(as.character(pupal_shape$area)) pupal_shape$area_mm <- (pupal_shape$area)*(1/(pupal_shape$scale_px/pupal_shape$scale_mm)^2) #Test if pupa area is different between the three species shapiro.test(pupal_shape$area_mm[pupal_shape$species=="suzukii_vincennes"]) #As it is not normal: test if Force depends on species kruskal.test(area_mm ~ species, data= pupal_shape) pairwise.wilcox.test(pupal_shape$area_mm , pupal_shape$species, p.adjust.method ="holm") #calculate median area pupa median(subset(pupal_shape$area_mm, pupal_shape$species == 'hydei_vincennes'), na.rm =T) median(subset(pupal_shape$area_mm, pupal_shape$species == 'simulans_vincennes'), na.rm =T) median(subset(pupal_shape$area_mm, pupal_shape$species == 'suzukii_vincennes'), na.rm =T) #calculate median area print / median area pupa median(subset(total_table_ok$Area_mm, total_table$Stock == "hydei_vincennes"),na.rm =T)/median(subset(pupal_shape$area_mm, pupal_shape$species == 'hydei_vincennes'), na.rm =T) median(subset(total_table_ok$Area_mm, total_table$Stock == "simulans_vincennes"),na.rm =T)/median(subset(pupal_shape$area_mm, pupal_shape$species == 'simulans_vincennes'), na.rm =T) median(subset(total_table_ok$Area_mm, total_table$Stock == "suzukii_vincennes"),na.rm =T)/median(subset(pupal_shape$area_mm, pupal_shape$species == 'suzukii_vincennes'), na.rm =T) #************************************************************************************************************************************************************************************* #************************************************************************************************************************************************************************************* #********************************************************************************************************************* ########################################################VINCENNES ASSAY######################################### rm(list=ls(all=T)) data_dif <- read.table("pupa_vincennes_final.txt", sep="\t", header=T) data_dif$Time <- factor(data_dif$Time, levels = c("t0", "t1", "t2", "t3", "t4", "t5" ,ordered = TRUE)) #install.packages("reshape2") library(reshape2) dat.m <- melt(data_dif,id.vars=c('bucket_ID','Orientation_C', 'Orientation_NC', 'Time'), measure.vars=c("count_NC", "count_C")) dat.m$bucket_ID <- as.character(dat.m$bucket_ID) #Boxplot with dots p3 <- ggplot(dat.m, aes(x=Time, y=value, fill= variable)) p3 <- p3 + geom_boxplot(width=0.5, outlier.shape = NA, position=position_dodge(0.9)) p3 <- p3 + geom_dotplot(alpha = 0.8, dotsize=0.35, binaxis='y', stackratio = 0.5,stackdir='center', position=position_dodge(0.9), colour= "black") p3 <- p3 + scale_y_continuous(name = "Number of pupae", breaks=c(0:10)) p3 <- p3 + theme_classic() p3 <- p3 + theme(axis.text=element_text(size=16), axis.title=element_text(size=16,face="bold")) p3 <- p3 + scale_fill_manual(values=c("white", "grey")) p3 #Statistics NC vs C #t1 wilcox.test(data_dif$count_C[data_dif$Time=="t1"], data_dif$count_NC[data_dif$Time=="t1"], paired=TRUE) #t2 wilcox.test(data_dif$count_C[data_dif$Time=="t2"], data_dif$count_NC[data_dif$Time=="t2"], paired=TRUE) #t3 wilcox.test(data_dif$count_C[data_dif$Time=="t3"], data_dif$count_NC[data_dif$Time=="t3"], paired=TRUE) #t4 wilcox.test(data_dif$count_C[data_dif$Time=="t4"], data_dif$count_NC[data_dif$Time=="t4"], paired=TRUE) #t5 wilcox.test(data_dif$count_C[data_dif$Time=="t5"], data_dif$count_NC[data_dif$Time=="t5"], paired=TRUE) ######################################################################################################## #************************************************************************************************************************************************************************************* #************************************************************************************************************************************************************************************* #********************************************************************************************************************* ########################################################PREDATION ASSAY IN THE LAB##################### rm(list=ls(all=T)) library(ggplot2) library(dplyr) library(data.table) library(gridExtra) ant <- read.table("pupa_predation_ant_nb.txt", sep="\t", header=T) data <- read.table("pupa_predation_results.txt", sep="\t", header=T) #************************************************** ##Calculate Duration of ants in contact with pupae data_m <- data #Add a column with unique ID for each pupa (colony and day of experiment + attached/detached) data_m$`variable condition` <- paste(data_m$variable, data_m$condition) #By summing all the time ant >0 #For attached ant_attached <- subset(ant, ant$condition == "attached") ant_attached <- ant_attached[,c(-1,-2)] duration_ant<- data.frame(colnames(ant_attached),"duration"=as.numeric(1)) for(i in 1:length(duration_ant[,1])) { duration_ant$duration[i] <- 5*(length(ant_attached[,i][ant_attached[i]>0]) - length(ant_attached[,i][is.na(ant_attached[i]==T)])) } duration_ant$colnames.ant_attached. <- paste(duration_ant$colnames.ant_attached.,"attached") #For detached ant_detached <- subset(ant, ant$condition == "detached") ant_detached <- ant_detached[,c(-1,-2)] duration_ant_2 <- data.frame(colnames(ant_detached),"duration"=as.numeric(1)) for(i in 1:length(duration_ant_2[,1])) { duration_ant_2$duration[i] <- 5*(length(ant_detached[,i][ant_detached[i]>0]) - length(ant_detached[,i][is.na(ant_detached[i]==T)])) } duration_ant_2$colnames.ant_detached. <- paste(duration_ant_2$colnames.ant_detached.,"detached") names(duration_ant)[1] <- "variable condition" names(duration_ant_2)[1] <- "variable condition" #bind duration_ant_f <- rbind(duration_ant, duration_ant_2) ##Merge both methods data_m <- merge(data_m, duration_ant_f, by="variable condition") #******************************************* ##Reformate table to do paired tests when needed data_m <- data.table(data_m) data.shape_m <- dcast.data.table(data_m, date + start_hour + end_hour + TRIAL_NB + colony_ID + variable + species ~ condition, value.var = c("strategies","condition","Time", "max_ant", "orientation","first_ant","duration","max_ant_at_first_done")) ############################################################################################# #********************************************************************************************** ##NOW USE data_m when long version of the table needed (Figure and some tests) and wide version of the #table when wide version of the table needed (paired test) #remove other data: rm(ant_attached, ant_detached, duration_ant, duration_ant_2, duration_ant_f) ############################################################################################# #********************************************************************************************** ###******************************************************** #ANT NB OVERTIME #**************************************************************** #Reshape ant table ant_t <- data.table(ant) ant.shape <- melt.data.table(ant_t, id.vars= c('Time', 'condition')) ant.shape$species <- NA #add species names to each experiment (Xa_b) with a: colony nb; b: trial nb levels(ant.shape$species) <- c(levels(ant.shape$species), "suzukii", "simulans") suz <- c('X1_1', "X3_1", 'X5_1', "X7_1", "X2_2", "X4_2", "X6_2", 'X1_3', "X3_3", 'X5_3', "X7_3", "X2_4", "X4_4", "X6_4",'X1_5', "X3_5", 'X5_5', "X7_5", "X2_6", "X4_6", "X6_6") sim <- c('X2_1', "X4_1", 'X6_1', "X1_2", "X3_2", "X5_2", "X7_2",'X2_3', "X4_3", 'X6_3', "X1_4", "X3_4", "X5_4", "X7_4",'X2_5', "X4_5", 'X6_5', "X1_6", "X3_6", "X5_6", "X7_6") for (i in 1:nrow(ant.shape)) { if (ant.shape$variable[i] %in% suz) { ant.shape$species[i] <- "suzukii" } else if (ant.shape$variable[i] %in% sim) { ant.shape$species[i] <- "simulans" } } ##Add bars for time to nest and time fully consumed p1 <- ggplot(data= ant.shape, aes(x= Time, y= value)) p1 <- p1 + geom_rect(aes(xmin = -10, xmax = 200, ymin = 0, ymax = 11, fill=species), alpha=0.8) #with lines p1 <- p1 + geom_vline(data= data_m, aes(xintercept= Time, color= condition, linetype = strategies), alpha=0.5) #With point #p1 <- p1 + geom_point(data= data, aes(x= Time, y=1, color= condition, shape = strategies)) p1 <- p1 + geom_line(aes(color= condition)) #Add Ants overtime p1 <- p1 + geom_hline(yintercept=0, alpha= 0.6) p1 <- p1 + facet_wrap(~ variable* species, ncol = 7 ) + theme_classic() + theme(legend.position = 'bottom') p1 <- p1 + scale_color_manual(values=c("darkred", "darkblue")) p1 <- p1 + scale_fill_manual(values= c("white", "lightgrey")) p1 <- p1 + theme(strip.text.x = element_text(margin = margin(b=2,t=2)), strip.text.y = element_text(margin = margin())) p1 ################################################################################### ###******************************************************** #ANT-PUPATION DURATION #**************************************************************** Fig_time <- data_m #Change Fig3 to time ants are found on pupae (duration or duration_dif) Fig_time$condition <- factor(Fig_time$condition , levels = c("detached", "attached", ordered = TRUE)) p2 <- ggplot(data= Fig_time, aes(x= condition, y= duration, fill= condition)) #p1 <- ggplot(data= Fig_time, aes(x= condition, y= duration_dif, fill= condition)) p2 <- p2 + geom_boxplot(outlier.shape = NA, width=0.8) p2 <- p2 + geom_point(alpha=0.4,shape = 21, colour = "black", size = 2, stroke = 1, position=position_jitter(0.02)) p2 <- p2 + theme_classic() + ylab("Ant-pupa duration (min)") + theme(plot.title = element_text(color="black", size=20, face = "bold"), axis.title = element_text(size = 20), axis.text = element_text(color="black", size=10)) p2 <- p2 + theme_classic() + theme(plot.title = element_text(color="black", size=20, face = "bold"), axis.title = element_text(size = 20), axis.text = element_text(color="black", size=12), strip.background = element_rect(fill = "white"), strip.text = element_text(size=15), legend.position =c(0.9, 0)) p2 <- p2 + scale_fill_manual(values=c("white", "grey")) p2 <- p2 + scale_y_continuous(breaks=seq(0,200,50)) p2 <- p2 + facet_wrap(~species) p2 #**************************************************************** #MAXIMUM ANT ON PUPAE over the experiment #**************************************************************** #Max ant on pupae with all the data p3 <- ggplot(data= Fig_time, aes(x= condition, y= max_ant, fill= condition)) p3 <- p3 + geom_boxplot(outlier.shape = NA, width=0.8) p3 <- p3 + geom_point(alpha=0.4,shape = 21, colour = "black", size = 2, stroke = 1, position=position_jitter(0.02)) p3 <- p3 + theme_classic() + ylab("Maximum ant on pupae") + theme(plot.title = element_text(color="black", size=20, face = "bold"), axis.title = element_text(size = 20), axis.text = element_text(color="black", size=10)) p3 <- p3 + theme_classic() + theme(plot.title = element_text(color="black", size=20, face = "bold"), axis.title = element_text(size = 20), axis.text = element_text(color="black", size=12), strip.background = element_rect(fill = "white"), strip.text = element_text(size=15), legend.position ='none') p3 <- p3 + scale_fill_manual(values=c("white", "grey")) p3 <- p3 + scale_y_continuous(limits = c(0, 11), breaks=seq(0,11,1)) p3 <- p3 + facet_wrap(~species) p3 #Calculate median duration of ant-pupa interaction for detached attached suz<- subset(Fig_time, Fig_time$species == "suzukii") sim<- subset(Fig_time, Fig_time$species == "simulans") dataMedian_suz <- summarise(group_by(suz, condition), MD = median(duration, na.rm = T)) dataMedian_sim <- summarise(group_by(sim, condition), MD = median(duration, na.rm = T)) #Calculate median max nb of ants for detached and attached pupae over the all experiment suz<- subset(Fig_time, Fig_time$species == "suzukii") sim<- subset(Fig_time, Fig_time$species == "simulans") dataMedian_suz <- summarise(group_by(suz, condition), MD = median(max_ant, na.rm = T)) dataMedian_sim <- summarise(group_by(sim, condition), MD = median(max_ant, na.rm = T)) #################### ##FIGURE #Fig3: 777 x 1000 px grid.arrange(p1,p2,p3,layout_matrix = rbind(c(1,1),c(1,1),c(2,3))) ############### ################################################ #########STATISTIC TESTS ######################Test pupa-ant duration detached vs attached (paired) #suzukii wilcox.test(data.shape_m$duration_detached[data.shape_m$species=="suzukii"], data.shape_m$duration_attached[data.shape_m$species=="suzukii"], paired=TRUE) #simulans wilcox.test(data.shape_m$duration_detached[data.shape_m$species=="simulans"], data.shape_m$duration_attached[data.shape_m$species=="simulans"], paired=TRUE) #Difference between attached and detached in both species ###Test pupa-ant duration suzukii vs simulans (not paired) #detached wilcox.test(data.shape_m$duration_detached[data.shape_m$species=="simulans"], data.shape_m$duration_detached[data.shape_m$species=="suzukii"]) #attached wilcox.test(data.shape_m$duration_attached[data.shape_m$species=="simulans"], data.shape_m$duration_attached[data.shape_m$species=="suzukii"]) #No difference between species ####################Test max ant attached vs detached (paired) #suzukii wilcox.test(data.shape_m$max_ant_detached[data.shape_m$species=="suzukii"], data.shape_m$max_ant_attached[data.shape_m$species=="suzukii"], paired=TRUE) #simulans wilcox.test(data.shape_m$max_ant_detached[data.shape_m$species=="simulans"], data.shape_m$max_ant_attached[data.shape_m$species=="simulans"], paired=TRUE) #Difference between attached and detached in both species ###Test max ant suzukii vs simulans (not paired) #detached wilcox.test(data.shape_m$max_ant_detached[data.shape_m$species=="simulans"], data.shape_m$max_ant_detached[data.shape_m$species=="suzukii"]) #attached wilcox.test(data.shape_m$max_ant_attached[data.shape_m$species=="simulans"], data.shape_m$max_ant_attached[data.shape_m$species=="suzukii"]) #No difference between species ##################Test if difference between max_ant_at_first_done detached vs attached #suzukii wilcox.test(data.shape_m$max_ant_at_first_done_detached[data.shape_m$species=="suzukii"], data.shape_m$max_ant_at_first_done_attached[data.shape_m$species=="suzukii"], paired=TRUE) #simulans wilcox.test(data.shape_m$max_ant_at_first_done_detached[data.shape_m$species=="simulans"], data.shape_m$max_ant_at_first_done_attached[data.shape_m$species=="simulans"], paired=TRUE) #Difference between attached and detached in both species #detached wilcox.test(data.shape_m$max_ant_at_first_done_detached[data.shape_m$species=="suzukii"], data.shape_m$max_ant_at_first_done_detached[data.shape_m$species=="simulans"]) #attached wilcox.test(data.shape_m$max_ant_at_first_done_attached[data.shape_m$species=="suzukii"], data.shape_m$max_ant_at_first_done_attached[data.shape_m$species=="simulans"]) #Difference between attached and detached in both species #######################Test first discovery time detached vs attached (paired) #suzukii wilcox.test(data.shape_m$first_ant_detached[data.shape_m$species=="suzukii"], data.shape_m$first_ant_attached[data.shape_m$species=="suzukii"], paired=TRUE) #simulans wilcox.test(data.shape_m$first_ant_detached[data.shape_m$species=="simulans"], data.shape_m$first_ant_attached[data.shape_m$species=="simulans"], paired=TRUE) #Difference between attached and detached only in suzukii #Test first discovery time suzukii vs simulans (not paired) #detached wilcox.test(data.shape_m$first_ant_detached[data.shape_m$species=="simulans"], data.shape_m$first_ant_detached[data.shape_m$species=="suzukii"]) #attached wilcox.test(data.shape_m$first_ant_attached[data.shape_m$species=="simulans"], data.shape_m$first_ant_attached[data.shape_m$species=="suzukii"]) #No difference between species ##########################Test Time to bring detached pupae to the nest suzukii VS simulans wilcox.test(data.shape_m$Time_detached[data.shape_m$strategies_detached == "to_nest" & data.shape_m$species=="simulans"], data.shape_m$Time_detached[data.shape_m$strategies_detached == "to_nest" & data.shape_m$species=="suzukii"]) #No difference between species ###################################################################################### ###METHODS PART ##################ADHESION PART rm(list=ls(all=T)) data <- read.table("pupa_adhesion_print.txt", sep="\t", header=T) #Stat on cuticle broke and not detached samples suz <- subset(data, data$Stock == "suzukii_vincennes") sim <- subset(data, data$Stock == "simulans_vincennes") hydei <- subset(data, data$Stock == "hydei_vincennes") length(subset(hydei, hydei$Comment_on_this_sample == "not_detached")[,1]) length(subset(hydei, hydei$Comment_on_this_sample == "cuticle_broke")[,1])