123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152115311541155115611571158115911601161116211631164116511661167116811691170117111721173117411751176117711781179118011811182118311841185118611871188118911901191119211931194119511961197119811991200120112021203120412051206120712081209121012111212121312141215121612171218121912201221122212231224122512261227122812291230123112321233123412351236123712381239124012411242124312441245124612471248124912501251125212531254125512561257125812591260126112621263126412651266126712681269127012711272127312741275127612771278127912801281128212831284128512861287128812891290129112921293129412951296129712981299130013011302130313041305130613071308130913101311131213131314131513161317131813191320132113221323132413251326132713281329133013311332133313341335133613371338133913401341134213431344134513461347134813491350135113521353135413551356135713581359136013611362136313641365136613671368136913701371137213731374137513761377137813791380138113821383138413851386138713881389139013911392139313941395139613971398139914001401140214031404140514061407140814091410141114121413141414151416141714181419142014211422142314241425142614271428142914301431143214331434143514361437143814391440144114421443144414451446144714481449145014511452145314541455145614571458145914601461146214631464146514661467146814691470147114721473147414751476147714781479148014811482148314841485148614871488148914901491149214931494149514961497149814991500150115021503150415051506150715081509151015111512151315141515151615171518151915201521152215231524152515261527152815291530153115321533153415351536153715381539154015411542154315441545154615471548154915501551155215531554155515561557155815591560156115621563156415651566156715681569157015711572157315741575157615771578157915801581158215831584158515861587158815891590159115921593159415951596159715981599160016011602160316041605160616071608160916101611161216131614161516161617161816191620162116221623162416251626162716281629163016311632163316341635163616371638163916401641164216431644164516461647164816491650165116521653165416551656165716581659166016611662166316641665166616671668166916701671167216731674167516761677167816791680168116821683168416851686168716881689169016911692169316941695169616971698169917001701170217031704170517061707170817091710171117121713171417151716171717181719172017211722172317241725172617271728172917301731173217331734173517361737173817391740174117421743174417451746 |
- ##' Data Analysis of Electroantennogram Data from two different Stonefly Ecotypes
- ##' (full-winged and wing-reduced)
- ##'
- ##' Main R routine to reproduce the Fig 1, Fig 2, and Fig S1 from the manuscript
- ##' entitled 'Rapid evolution of olfactory degradation in recently flightless
- ##' insects'.
- ##'
- ##
- # load R packages
- library(dplyr)
- library(tidyr)
- library(ggplot2)
- library(tibble)
- library(ggbeeswarm)
- library(scales)
- library(arm)
- library(blmeco)
- library(lemon)
- library(ggpubr)
- library(grid)
- library(gridExtra)
- library(magick)
- library(cowplot)
- library(grImport)
- library(pdftools)
- library(multitaper)
- # clean up environment
- rm(list=ls())
- # set working directory
- setwd(dirname(rstudioapi::getActiveDocumentContext()$path))
- # load function to calculate various response dynamic measures
- source("./functions/responseDynamics.R")
- source('./functions/utils.R')
- fig_path <- "../Figures"
- ### colors for stoneflies ###
- fw_l <- "#787878" ## full-winged: grey
- wr_l <- "#00CC00" ## wing-reduced: green
- ## load the sub plots for Fig. 1
- Fig1A <- magick::image_read_pdf('../Pictures/Fig1A.pdf')
- Fig1A <- ggplot() +
- theme_bw() +
- theme(panel.border = element_blank()) +
- draw_image(Fig1A)
- Fig1B <- magick::image_read_pdf('../Pictures/Fig1B.pdf')
- Fig1B <- ggplot() +
- theme_bw() +
- theme(panel.border = element_blank()) +
- draw_image(Fig1B)
- # number of simulation draws from posterior distribution
- nsim <- 10000
- # aspect ratio (as y/x) for trace plots
- ar <- 2/3
- ## Read data set
- data2 <- read.csv('../Datasets/Dataset1.csv', stringsAsFactors = FALSE)
- # create median filtered mean traces (only for non-10Hz stimulations)
- med <- as.data.frame(t(apply(X = data2[data2$PulseDuration != 50, paste0("t", 1:5000)],
- MARGIN = 1,
- FUN = runmed, k = 11)))
- data2[data2$PulseDuration != 50, paste0("t", 1:5000)] <- med
- ### automatically exclude alive Antennae that are not responding to 2-HeptanonTest ######
- heptTest <- data2 %>%
- filter(Odor=="2HeptanoneTest", AntennaState=="alive")
- par(mfrow=c(4,4))
- analyze <- vector("logical", 0)
- heptTest$keep <- NA
- for(i in 1:nrow(heptTest)) {
- # in which ms is the signal higher than 2 x sd during 2 to 5 sec win of trace
- num <- which(heptTest[i, paste0("t", 250:1200)] > 2 * sd(heptTest[i, paste0("t", 2000:5000)]))
- result <- rle(diff(num))
- ## above threshold for at least 50 ms in a row?
- an <- any(result$lengths >= 50 & result$values==1)
- heptTest[i, "keep"] <- an
- analyze <- c(analyze, an)
- }
- ID_keep <- unique(heptTest$ID_new[heptTest$keep])
- ID_keep <- c(ID_keep, unique(data2[data2$Group=="PID", "ID_new"]))
- data2 <- data2[data2$ID_new %in% ID_keep, ]
- data2 <- droplevels(data2)
- ## Correct for amplification
- ## Note: 1 unit for different groups has different voltage ranges
- ## PID - 30 V
- ## Bee - 10 mV
- ## Fenestrata - 1 mV
- multiply_voltage <- 1/data2$Amplification
- data2[, paste0("t", 1:5000)] <- data2[, paste0("t", 1:5000)] * multiply_voltage # in Volt
- data2[data2$Group != "PID", paste0("t", 1:5000)] <- data2[data2$Group != "PID", paste0("t", 1:5000)] * 1000 # in millivolt
- ## Convert some columns to a factor
- data2$ID_new <- as.factor(data2$ID_new)
- data2$Winged <- factor(data2$Winged, levels = c("full-winged", "wing-reduced", "none"))
- data2$Location <- as.character(data2$Location)
- data2$Location[data2$Group=="PID"] <- "PID"
- data2$Location[data2$Species=="Honeybee"] <- "Bee"
- # reorder the factor levels of column Odor and adjust the names
- data2$Odor <- factor(data2$Odor, levels = c("Blank", "2Heptanone", "1Octanol",
- "PropionicAcid", "2Butanone", "Stonefly", "2HeptanoneTest"),
- labels = c("Blank","2-heptanone", "1-octanol",
- "Propionic acid", "2-butanone", "Stonefly", "2-heptanone (2nd)"))
- data2$Group <- factor(data2$Group, levels = c("Burns",
- "BlackJack", "Whisky", "SixMile",
- "Lug", "Bee", "PID"),
- labels = c("Mt Burns", "Black Jacks", "Whiskey", "Six Mile", "Lug", "Bee", "PID"))
- data2 <- data2 %>%
- filter(AntennaState == "alive" | Group == "PID")
- data2$NewGroup <- paste(data2$Group, data2$Winged)
- data2$NewGroup <- as.factor(data2$NewGroup)
- data2 <- droplevels(data2)
- new_colnames <- c("Timeto10Max",
- "Timeto90Max",
- "TimetoMax",
- "TimeFromMaxto90",
- "TimeFromMaxto10",
- "Time10to90Max",
- "Time90to10Max",
- "Timeto10Min",
- "Timeto90Min",
- "TimetoMin",
- "TimeFromMinto90",
- "TimeFromMinto10",
- "Time10to90Min",
- "Time90to10Min",
- "MeanMax",
- "MeanMin")
- y1 <- as.data.frame(t(apply(data2[,paste0("t", 201:5000)], 1, responseDynamics, from=0.1, to=0.9)))
- colnames(y1) <- new_colnames
- data2 <- cbind(as.data.frame(data2), y1)
- data2 <- data2 %>%
- mutate(Timeto10 = ifelse(Odor == "Propionic acid" |
- Odor == "Stonefly", Timeto10Min, Timeto10Max),
- Time10to90 = ifelse(Odor == "Propionic acid" |
- Odor == "Stonefly", Time10to90Min, Time10to90Max),
- Time90to10 = ifelse(Odor == "Propionic acid" |
- Odor == "Stonefly", Time10to90Min, Time90to10Max))
- data2$max <- apply(data2[,paste0("t", 251:1200)], 1, max)
- data2$min <- apply(data2[,paste0("t", 251:1200)], 1, min)
- data2 <- data2 %>%
- mutate(mima = ifelse(Odor == "Propionic acid" |
- Odor == "Stonefly", min, max))
- data2 <- data2 %>%
- mutate(Mean = ifelse(Odor == "Propionic acid" |
- Odor == "Stonefly", MeanMin, MeanMax))
- data2 <- data2 %>%
- mutate(TimeBackTo10 = ifelse(Odor == "Propionic acid" |
- Odor == "Stonefly", TimeFromMinto10, TimeFromMaxto10))
- data2 %>%
- filter(Odor=="2-heptanone" | Odor=="2-heptanone (2nd)", PulseDuration==300) %>%
- group_by(Winged, ID_new, AntennaState) %>%
- summarise(ratio = mima[Odor=="2-heptanone (2nd)"] / mima[Odor=="2-heptanone"]) %>%
- group_by(Winged, AntennaState) %>%
- summarise(mean=mean(ratio), sd=sd(ratio), median=median(ratio), n=n())
- #######################################
- ### How big is the sample size for the different groups?
- data2 %>%
- group_by(Location, Group, Winged, ID_new) %>%
- summarise(n=n()) %>%
- group_by(Location, Group, Winged) %>%
- summarise(n=n())
- tst <- data2
- tst <- droplevels(tst)
- table(tst$Config, tst$Group, tst$Winged)
- locations <- c("Mt Burns", "Black Jacks", "Whiskey", "Six Mile", "Lug")
- odors <- c("2-heptanone", "1-octanol", "2-butanone", "2-heptanone (2nd)")
- ############## Statistics #############
- ## Bayesian Statistics ###
- creeks <- locations
- plot_bayes_mima <- function(data, groups) {
-
- ################################################################################
-
- data1 <- data2
-
- odors <- c("2-heptanone", "1-octanol", "2-butanone")
-
- df <- as.data.frame(matrix(NA, ncol=7, nrow = 0))
- newDf <- as.data.frame(matrix(NA, ncol=3, nrow = 0))
- newDf_all <- as.data.frame(matrix(NA, ncol=3, nrow = 0))
-
- for (odor in odors) {
- if (odor == '2-heptanone') {
- data <- data1 %>% filter(Odor %in% c('2-heptanone (2nd)', '2-heptanone'))
- } else {
- data <- data1 %>% filter(Odor == odor)
- }
-
- data <- data %>% filter(Group %in% groups, PulseDuration != 50)
-
- data$PulseDuration <- as.character(data$PulseDuration)
-
- data$PulseDuration <- ifelse(data$PulseDuration == 300 & data$Odor == '2-heptanone (2nd)',
- '300\n (2nd)', data$PulseDuration)
- data$Odor <- as.character(data$Odor)
-
- data$Odor <- ifelse(data$Odor == '2-heptanone (2nd)', '2-heptanone', data$Odor)
-
- data$Location <- factor(data$Location)
- data$PulseDuration <- factor(data$PulseDuration, levels = c('15', '30', '150', '300', '300\n (2nd)'))
- data <- droplevels(data)
-
-
-
-
- for (cc in seq_along(groups)) {
- subdata <- data %>% filter(Group == groups[cc])
- subdata$PulseDuration <- as.factor(subdata$PulseDuration)
- subdata <- droplevels(subdata)
- min_mean <- min(subdata$Mean)
-
- if (min_mean < 0) {
- subdata$adj_Mean <- subdata$Mean + abs(min_mean) + 0.01
- } else {
- subdata$adj_Mean <- subdata$Mean
- }
-
- subdata$log2 <- log2(subdata$adj_Mean)
-
- if (nrow(subdata) == 0)
- next
-
- if (length(unique(subdata$PulseDuration)) == 1) {
-
- m1 <- lm(log2 ~ Winged, data=subdata)
-
- } else {
-
- if (unique(subdata$Winged == 'none')) {
-
- m1 <- lmer(log2 ~ PulseDuration + (1 | ID_new), data=subdata, REML = FALSE)
-
- } else {
-
- m1 <- lmer(log2 ~ Winged * PulseDuration + (1 | ID_new), data=subdata, REML = FALSE)
-
- }
-
-
- }
-
-
- par(mfrow=c(2,2))
-
- if (length(unique(subdata$PulseDuration)) != 1) {
- scatter.smooth(fitted(m1), resid(m1))
- abline(h=0, lty=2)
- title("Tukey-Anscombe Plot") # residuals vs. fitted
-
- qqnorm(resid(m1), main="normal QQ-plot, residuals") # qq of residuals
- qqline(resid(m1))
- scatter.smooth(fitted(m1), sqrt(abs(resid(m1)))) # res. var vs. fitted
-
- qqnorm(ranef(m1)$ID_new[,1], main="normal QQ-plot, random effects")
- qqline(ranef(m1)$ID_new[,1]) # qq of random effects
- } else {
- plot(m1)
- }
-
-
- summary(m1)
-
-
-
-
-
- if (length(unique(subdata$PulseDuration)) == 1) {
-
- tmp <- unique(paste( subdata$PulseDuration, subdata$Winged, subdata$Odor, sep = '_'))
-
- tmp <- as.data.frame(matrix(unlist(strsplit(tmp, split = '_')), ncol=3, byrow = TRUE))
-
- tmp <- tmp[order(as.numeric(as.character(tmp$V1)), tmp$V2, tmp$V3), ]
- colnames(tmp) <- c('PulseDuration', 'Winged', 'Odor')
-
- tmp$PulseDuration <- factor(tmp$PulseDuration)
- tmp$PulseDuration <- factor(tmp$PulseDuration,
- levels = levels(tmp$PulseDuration)[order(as.numeric(as.character(levels(tmp$PulseDuration))))])
-
- tmp$Winged <- factor(tmp$Winged, levels = c('full-winged', 'wing-reduced'))
- tmp$Odor <- factor(tmp$Odor, levels = c('2-heptanone', '1-octanol', '2-butanone'))
-
- newdat <- tmp
-
- bsim <- sim(m1, n.sim=nsim)
-
- fitmat <- matrix(ncol=nsim, nrow=nrow(newdat))
- Xmat <- model.matrix( ~ Winged, data=newdat)
-
-
- for(i in 1:nsim) fitmat[,i] <- Xmat %*% bsim@coef[i,] #fitted values
-
-
- tmp <- t(apply(fitmat, 1, quantile, prob=c(0.025, 0.5, 0.975)))
-
- if (min_mean < 0) {
-
- tmp <- 2^(tmp) - abs(min_mean) - 0.01
-
- } else {
-
- tmp <- 2^(tmp)
-
- }
-
-
-
-
- tst <- cbind(newdat, tmp, 'Group' = groups[cc])
-
- df <- rbind(df, tst)
-
- for (oo in levels(newdat$Odor)) {
-
- datt <- newdat[newdat$Odor == oo, ]
-
- for (pp in unique(datt$PulseDuration)) {
-
- dat <- fitmat[which(newdat$PulseDuration == pp & newdat$Odor == oo), ]
- star2 <- sum(dat[1, ] > dat[2, ]) / nsim
-
- newDf <- rbind(newDf,
- data.frame('Group' = groups[cc],
- 'Odor' = oo,
- 'Diff' = star2,
- 'PulseDuration' = pp))
-
- }
-
-
- }
-
-
- } else {
-
- if (unique(subdata$Winged == 'none')) {
-
- #######################
-
- tmp <- unique(paste(subdata$PulseDuration, subdata$Winged, subdata$Odor, sep = '_'))
-
- tmp <- as.data.frame(matrix(unlist(strsplit(tmp, split = '_')), ncol=3, byrow = TRUE))
-
- tmp <- tmp[order(as.numeric(as.character(tmp$V1)), tmp$V2, tmp$V3), ]
- colnames(tmp) <- c('PulseDuration', 'Winged', 'Odor')
-
- tmp$PulseDuration <- factor(tmp$PulseDuration)
- tmp$PulseDuration <- factor(tmp$PulseDuration,
- levels = levels(tmp$PulseDuration)[order(as.numeric(as.character(levels(tmp$PulseDuration))))])
- tmp$Winged <- factor(tmp$Winged, levels = c('full-winged', 'wing-reduced', 'none'))
- tmp$Odor <- factor(tmp$Odor, levels = c('2-heptanone', '1-octanol', '2-butanone'))
-
-
- newdat <- tmp
-
-
- bsim <- sim(m1, n.sim=nsim)
-
- fitmat <- matrix(ncol=nsim, nrow=nrow(newdat))
- Xmat <- model.matrix( ~ PulseDuration, data=newdat)
-
-
- for(i in 1:nsim) fitmat[,i] <- Xmat %*% bsim@fixef[i,] #fitted values
-
- names(Xmat)
- length(newdat$Odor)
-
- tmp <- t(apply(fitmat, 1, quantile, prob=c(0.025, 0.5, 0.975)))
-
- if (min_mean < 0) {
-
- tmp <- 2^(tmp) - abs(min_mean) - 0.01
-
- } else {
-
- tmp <- 2^(tmp)
-
- }
-
- tst <- cbind(newdat, tmp, 'Group' = groups[cc])
-
- df <- rbind(df, tst)
-
-
- for (oo in levels(newdat$Odor)) {
-
- datt <- newdat[newdat$Odor == oo, ]
-
- for (pp in unique(datt$PulseDuration)) {
-
- dat <- fitmat[which(newdat$PulseDuration == pp & newdat$Odor == oo), ]
- newDf <- rbind(newDf,
- data.frame('Group' = groups[cc],
- 'Odor' = oo,
- 'Diff' = 0.5,
- 'PulseDuration' = pp))
-
- }
-
-
- }
-
- #######################
-
- } else {
-
- tmp <- unique(paste( subdata$PulseDuration, subdata$Winged, subdata$Odor, sep = '_'))
-
- tmp <- as.data.frame(matrix(unlist(strsplit(tmp, split = '_')), ncol=3, byrow = TRUE))
-
- tmp <- tmp[order(as.numeric(as.character(tmp$V1)), tmp$V2, tmp$V3), ]
- colnames(tmp) <- c('PulseDuration', 'Winged', 'Odor')
-
- tmp$PulseDuration <- factor(tmp$PulseDuration)
- tmp$PulseDuration <- factor(tmp$PulseDuration,
- levels = levels(tmp$PulseDuration)[order(as.numeric(as.character(levels(tmp$PulseDuration))))])
- tmp$Winged <- factor(tmp$Winged, levels = c('full-winged', 'wing-reduced'))
- tmp$Odor <- factor(tmp$Odor, levels = c('2-heptanone', '1-octanol', '2-butanone'))
-
-
- newdat <- tmp
-
-
- bsim <- sim(m1, n.sim=nsim)
-
- fitmat <- matrix(ncol=nsim, nrow=nrow(newdat))
- Xmat <- model.matrix( ~ Winged * PulseDuration, data=newdat)
-
-
-
- for(i in 1:nsim) fitmat[,i] <- Xmat %*% bsim@fixef[i,] #fitted values
-
- tmp <- t(apply(fitmat, 1, quantile, prob=c(0.025, 0.5, 0.975)))
-
- if (min_mean < 0) {
-
- tmp <- 2^(tmp) - abs(min_mean) - 0.01
-
- } else {
-
- tmp <- 2^(tmp)
-
- }
-
- tst <- cbind(newdat, tmp, 'Group' = groups[cc])
-
- df <- rbind(df, tst)
-
-
- for (oo in levels(newdat$Odor)) {
-
- datt <- newdat[newdat$Odor == oo, ]
-
- for (pp in unique(datt$PulseDuration)) {
-
- dat <- fitmat[which(newdat$PulseDuration == pp & newdat$Odor == oo), ]
- star2 <- sum(dat[1, ] > dat[2, ]) / nsim
-
- newDf <- rbind(newDf,
- data.frame('Group' = groups[cc],
- 'Odor' = oo,
- 'Diff' = star2,
- 'PulseDuration' = pp))
-
- }
-
- }
-
- }
-
- }
-
- }
-
- }
-
-
- colnames(df)[4:6] <- c("lower", "fitted", "upper")
- newDf
-
- newDf$Star1 <- newDf$Diff >= 0.95 | newDf$Diff <= 0.05
- newDf$Star2 <- newDf$Diff >= 0.99 | newDf$Diff <= 0.01
- newDf$Star <- ''
- newDf$Star[newDf$Star1 == TRUE] <- '*'
- newDf$Star[newDf$Star2 == TRUE] <- '**'
-
- df$Group
- newDf$Winged <- 'star'
- col <- c(fw_l, wr_l, "blue")
-
- plotdata <- data1 %>%
- filter(Odor %in% c("2-heptanone", "2-heptanone (2nd)"),
- PulseDuration != 50)
- plotdata$PulseDuration <- as.character(plotdata$PulseDuration)
- plotdata$PulseDuration[plotdata$Odor == "2-heptanone (2nd)"] <- '300\n (2nd)'
-
- plotdata$Odor[plotdata$Odor == "2-heptanone (2nd)"] <- "2-heptanone"
-
- plotdata$PulseDuration <- factor(plotdata$PulseDuration, levels = c(15, 30, 150, 300, "300\n (2nd)"))
-
-
- plotdf <- df %>%
- filter(Odor %in% c("2-heptanone", "2-heptanone (2nd)"))
- plotdf$PulseDuration <- as.character(plotdf$PulseDuration)
- plotdf$PulseDuration[plotdf$Odor == "2-heptanone (2nd)"] <- '300\n (2nd)'
- plotdf$Odor[plotdf$Odor == "2-heptanone (2nd)"] <- "2-heptanone"
- plotdf$PulseDuration <- factor(plotdf$PulseDuration, levels = c(15, 30, 150, 300, "300\n (2nd)"))
-
-
- newDf_plot <- newDf %>%
- filter(Odor %in% c("2-heptanone", "2-heptanone (2nd)"))
- newDf_plot$PulseDuration <- as.character(newDf_plot$PulseDuration)
- newDf_plot$PulseDuration[newDf_plot$Odor == "2-heptanone (2nd)"] <- '300\n (2nd)'
- newDf_plot$Odor[newDf_plot$Odor == "2-heptanone (2nd)"] <- "2-heptanone"
- newDf_plot$PulseDuration <- factor(newDf_plot$PulseDuration, levels = c(15, 30, 150, 300, "300\n (2nd)"))
-
-
- fen <- ggplot(plotdata %>% filter(Group %in% creeks), aes(PulseDuration, Mean, color=Winged)) +
- facet_grid(Odor ~ Group) +
- geom_beeswarm(dodge.width = 0.5, pch = 1) +
- geom_errorbar(data=plotdf %>% filter(Group %in% creeks), aes(ymin=lower, ymax=upper, y=NULL),
- color="black", width=0.5, position = position_dodge2(padding = 1, width=0.1)) +
- stat_summary(data=newDf_plot %>% filter(Group %in% creeks),
- aes(x=as.factor(PulseDuration), y=0.15, label = Star),
- geom = 'text', fun = max, color="black",
- size = 8) +
- geom_beeswarm(data=plotdf %>% filter(Group %in% creeks), aes(y = fitted, group = Winged),
- pch = '_', dodge.width = 0.5, color = 'black', size=3) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- theme_bw() +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_line(colour = "black"),
- axis.title = element_blank(),
- strip.background = element_rect(colour="NA", fill="NA"),
- strip.text = element_blank(),
- panel.border = element_blank(),
- text = element_text(size=16),
- axis.text = element_text(colour="black"))
- fen
-
- bee <- ggplot(plotdata %>% filter(Group %in% 'Bee'), aes(PulseDuration, Mean, color=Winged)) +
- facet_grid(Odor ~ Group) +
- geom_beeswarm(dodge.width = 0.5, pch = 1) +
- geom_errorbar(data=plotdf %>% filter(Group %in% 'Bee'), aes(ymin=lower, ymax=upper, y=NULL), color="black",
- width=0, position = position_dodge2(padding = 1, width=0.1)) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = 'blue') +
- theme_bw() +
- geom_beeswarm(data=plotdf %>% filter(Group %in% 'Bee'), aes(y=fitted, group=Winged),
- pch = '_', color = 'black', groupOnX = TRUE, size=3) +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_line(colour = "black"),
- axis.title=element_blank(),
- strip.background = element_rect(colour="NA", fill="NA"),
- strip.text = element_blank(),
- panel.border = element_blank(),
- text = element_text(size=16),
- axis.text=element_text(colour="black"))
- bee
-
- pid <- ggplot(plotdata %>% filter(Group %in% 'PID'), aes(PulseDuration, Mean, color=Winged)) +
- facet_grid(Odor ~ Group) +
- geom_beeswarm(dodge.width = 0.5, pch = 1) +
- geom_errorbar(data=plotdf %>% filter(Group %in% 'PID'), aes(ymin=lower, ymax=upper, y=NULL),
- color="black", width=0, position = position_dodge2(padding = 1, width=0.1)) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = 'blue') +
- theme_bw() +
- geom_beeswarm(data=plotdf %>% filter(Group %in% 'PID'), aes(y=fitted, group=Winged),
- pch = '_', dodge.width = 0.5, color = 'black', groupOnX = TRUE, size=3) +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_line(colour = "black"),
- axis.title=element_blank(),
- plot.tag = element_text(size = 14, hjust = 0.5, vjust = 0.5),
- plot.tag.position = c(0.02, 1 - 0.03),
- strip.background = element_rect(colour="NA", fill="NA"),
- strip.text = element_blank(),
- panel.border = element_blank(),
- text = element_text(size=16),
- axis.text=element_text(colour="black"))
- pid
-
- l <- list('Diffs' = newDf, 'Mean_crI' = df, 'fen' = fen, 'bee' = bee, 'pid' = pid)
- return(l)
- }
- ### Figure S1 #####
- Fly_mean <- plot_bayes_mima(data = data2, groups = c(creeks, 'Bee', 'PID'))
- FigS1B <- egg::ggarrange(Fly_mean$fen, Fly_mean$bee, Fly_mean$pid, clip = FALSE,
- labels = c('', '', expression(x ~ 10^3)), label.args = list(vjust=0),
- nrow = 1, ncol=3, widths = c(5, 1, 1), heights = c(6), debug = FALSE)
- y.grob <- textGrob(paste("Response strength (mV)"),
- gp=gpar(fontsize=18), rot=90)
- x.grob <- textGrob("Pulse duration (ms)",
- gp=gpar(fontsize=18))
- #add to plot
- FigS1B <- grid.arrange(arrangeGrob(FigS1B, left = y.grob, bottom = x.grob))
- plot_winged_all <- function(data, odors, creeks, Diffs) {
-
- data <- data2
- data <- data %>%
- filter(Group %in% c(creeks, "Bee", "PID")) %>%
- filter(PulseDuration == 300) %>%
- filter(Odor %in% odors)
-
- data <- droplevels(data)
-
- data$Odor <- factor(data$Odor, levels = levels(data$Odor),
- labels = c("2-heptanone", "1-octanol", "2-butanone", "2-heptanone\n (2nd)"))
-
- # calculate the mean over all odor-pulseduration configurations
- d_mean <- data %>%
- group_by(Group, Winged, Odor, PulseDuration) %>%
- summarise_at(.funs = mean, .vars=paste0("t", 1:5000))
-
-
- # calculate the standard error over all odor-pulseduration configurations
- d_se <- data %>%
- group_by(Group, Winged, Odor, PulseDuration) %>%
- summarise_at(.funs = se, .vars=paste0("t", 1:5000))
-
- # calculate sample size (n) over all odor-pulseduration configurations
- d_n <- data %>%
- group_by(Group, Winged, Odor, PulseDuration) %>%
- summarise(N=n()) %>%
- group_by(Group, Winged) %>%
- summarise(N=unique(N))
-
- d_max <- d_mean %>%
- group_by(Group, Odor, PulseDuration) %>%
- summarise(Max=max(t400))
-
-
- # create d_all and add the sample size
- d_all <- d_mean
- data_se <- gather(d_se, time, se, t1:t5000, factor_key=TRUE)
-
-
- d_all$PulseDuration <- factor(d_all$PulseDuration, levels = c(300))
- d_all$PulseDur <- factor(d_all$PulseDuration, levels = c(300),
- labels = c("300 ms"))
-
-
- data_long <- gather(d_all, time, Voltage, t1:t5000, factor_key=TRUE)
- data_long$time <- (as.numeric(gsub("t", "", data_long$time)) - 200)/1000
- data_long$SE <- data_se$se
-
-
- data_long$high <- data_long$Voltage + data_long$SE
- data_long$low <- data_long$Voltage - data_long$SE
-
-
- tmp <- expand.grid(creeks, c("FW", "WR"), stringsAsFactors = FALSE)
- tmp <- tmp[order(tmp$Var1), ]
- tmp <- paste0(tmp$Var2, ", ", tmp$Var1)
-
- tmp <- c(tmp, "Bee", "PID")
-
- data_long$NewGroup <- factor(paste(data_long$Winged, data_long$Group),
- labels = tmp)
-
- col <- c(fw_l, wr_l, "blue")
-
-
- print(d_n)
-
- geom.text.size <- 5.2
- data_long$Winged <- factor(data_long$Winged,
- levels = levels(data_long$Winged),
- labels = c("Full-winged", "Wing-reduced", "none"))
-
-
- data_segm1 <- data.frame(y=0, yend=0.1, x=Inf, xend=Inf,
- PulseDuration= c(300),
- Group="Lug",
- Odor=c("2-heptanone\n (2nd)"))
- label1 <- data.frame(x=Inf, y=0.05,
- text=c("0.1 mV"),
- PulseDuration= c(300),
- Group="Lug",
- Odor=c("2-heptanone\n (2nd)"))
-
- fw_text <- data.frame(x = 1, y = 0.07,
- text=c("full-winged"),
- PulseDuration= c(300),
- Group="Mt Burns",
- Odor=c("2-heptanone"))
- wr_text <- data.frame(x = 1, y = 0.1,
- text=c("wing-reduced"),
- PulseDuration= c(300),
- Group="Mt Burns",
- Odor=c("2-heptanone"))
-
- Diffs$Odor <- as.character(Diffs$Odor)
- Diffs$Odor[Diffs$PulseDuration == '300\n (2nd)'] <- '2-heptanone\n (2nd)'
- Diffs$Odor <- factor(Diffs$Odor)
- Diffs$Odor <- factor(Diffs$Odor, levels = c('2-heptanone', '1-octanol', '2-butanone', '2-heptanone\n (2nd)'))
-
- Diffs$PulseDuration[Diffs$PulseDuration == '300\n (2nd)'] <- '300'
-
- Diffs <- droplevels(Diffs)
-
- Diffs <- merge(Diffs, d_max, by = c('Group', 'Odor', 'PulseDuration'))
-
- ## plotting mean and se
- Fen <- ggplot(data_long %>% filter(Group %in% creeks), aes(x=time, y=Voltage)) +
- facet_grid(Odor ~ Group, switch = "y") +
- geom_rect(mapping=aes(xmin=0, xmax=as.numeric(as.character(PulseDuration))/1000,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_ribbon(aes(ymin=Voltage + SE, ymax=Voltage - SE, fill=Winged),
- color=NA, alpha=0.3) +
- geom_line(aes(colour = Winged)) +
- theme_bw() +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_blank(),
- axis.title = element_blank(),
- plot.margin = unit(c(5.5, 5.5, 5.5, 5.5), "points"),
- strip.background = element_rect(colour="NA", fill="NA"),
- panel.border = element_blank(),
- text = element_text(size=18),
- axis.ticks = element_blank(),
- aspect.ratio = ar * 4,
- strip.text.y.left = element_text(angle = 0, size = 16),
- strip.text.x = element_text(size = 16),
- axis.text=element_blank()) +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- geom_text(data = wr_text, aes(x=x, y=y, label=text), size=4, hjust=0, color=wr_l) +
- geom_text(data = fw_text, aes(x=x, y=y, label=text), size=4, hjust=0, color=fw_l) +
- geom_text(data = Diffs, aes(x=0.4, y=Max + 0.025, label=Star), size=geom.text.size, hjust=0.5) +
- coord_cartesian(expand = TRUE,
- clip = 'off')
- Fen
-
- data_segm1 <- data.frame(y=0, yend=1, x=Inf, xend=Inf,
- PulseDuration= c(300),
- Group="Bee",
- Odor=c("2-heptanone\n (2nd)"))
- label1 <- data.frame(x=Inf, y=0.5,
- text=c("1 mV"),
- PulseDuration= c(300),
- Group="Bee",
- Odor=c("2-heptanone\n (2nd)"))
- seg <- data.frame(x=4.5, y = -0.35,
- text=c(""),
- PulseDuration= c(300),
- Group="Bee",
- Odor=c("2-heptanone\n (2nd)"))
-
-
- B <- ggplot(data_long %>% filter(Group %in% "Bee"), aes(x=time, y=Voltage)) +
- facet_grid(Odor ~ Group, switch = "y") +
- geom_rect(mapping=aes(xmin=0, xmax=as.numeric(as.character(PulseDuration))/1000,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = "Blue") +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_ribbon(aes(ymin=Voltage + SE, ymax=Voltage - SE, fill=Winged),
- color=NA, alpha=0.3) +
- geom_line(aes(colour = Winged)) +
- theme_bw() +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_blank(),
- axis.title = element_blank(),
- plot.margin = unit(c(5.5, 5.5, 5.5, 5.5), "points"),
- strip.background = element_rect(colour="NA", fill="NA"),
- panel.border = element_blank(),
- text = element_text(size=18),
- axis.ticks = element_blank(),
- aspect.ratio = ar * 4,
- strip.text.y.left = element_blank(),
- strip.text.x = element_text(size = 16),
- axis.text=element_blank()) +
- geom_blank(data = seg, mapping = aes(x=x, y=y), inherit.aes = FALSE) +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- coord_cartesian(expand = TRUE,
- clip = 'off')
-
- data_segm1 <- data.frame(y=0, yend=c(3, 0.05, 3), x=Inf, xend=Inf,
- PulseDuration= c(300),
- Group="PID",
- Odor=c("2-heptanone\n (2nd)", "1-octanol", "2-heptanone"))
- label1 <- data.frame(x=Inf, y=c(1.5, 0.04, 1.5),
- text=c("3 V", "0.05 V", "3 V"),
- PulseDuration= c(300),
- Group="PID",
- Odor=c("2-heptanone\n (2nd)", "1-octanol","2-heptanone"))
-
- data_segm2 <- data.frame(x=3.5, xend=4.5, y=-0.3, yend=-0.3,
- PulseDuration= c(300),
- Group="PID",
- Odor=c("2-heptanone\n (2nd)"))
-
- label2 <- data.frame(x=4, y=-0.2,
- text=c("1 s"),
- PulseDuration= c(300),
- Group="PID",
- Odor=c("2-heptanone\n (2nd)"))
-
- seg <- data.frame(x=4.5, y = -0.35,
- text=c(""),
- PulseDuration= c(300),
- Group="PID",
- Odor=c("2-butanone"))
-
- P <- ggplot(data_long %>% filter(Group %in% "PID"), aes(x=time, y=Voltage)) +
- facet_grid(Odor ~ Group, switch = "y", scales = 'free_y') +
- geom_rect(mapping=aes(xmin=0, xmax=as.numeric(as.character(PulseDuration))/1000,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = "Blue") +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_ribbon(aes(ymin=Voltage + SE, ymax=Voltage - SE, fill=Winged),
- color=NA, alpha=0.3) +
- geom_line(aes(colour = Winged)) +
- theme_bw() +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_blank(),
- axis.title = element_blank(),
- plot.margin = unit(c(5.5, 5.5, 5.5, 5.5), "points"),
- strip.background = element_rect(colour="NA", fill="NA"),
- panel.border = element_blank(),
- text = element_text(size=18),
- axis.ticks = element_blank(),
- aspect.ratio = ar * 4,
- strip.text.y.left = element_blank(),
- strip.text.x = element_text(size = 16),
- axis.text=element_blank()) +
- geom_blank(data = seg, mapping = aes(x=x, y=y), inherit.aes = FALSE) +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- geom_segment(data=data_segm2,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label2, aes(x=x, y=y, label=text), size=geom.text.size, vjust=1.3) +
- coord_cartesian(expand = TRUE,
- clip = 'off')
- P
-
- cowplot::plot_grid(Fen, B, P, nrow = 1, align = "h")
- list1 <- list("traceW"=Fen,
- "traceB"=B,
- "traceP"=P)
- }
- traces300 <- plot_winged_all(data = data2,
- odors = odors,
- creeks = locations,
- Diffs = Fly_mean$Diffs %>% filter(Group %in% creeks,
- PulseDuration %in% c('300', '300\n (2nd)')))
- Fig1C <- egg::ggarrange(traces300$traceW, traces300$traceB, traces300$traceP,
- heights = 5, nrow = 1, ncol=3, widths = c(5, 1, 1), debug = FALSE)
- ########################
- data <- data2
- data <- data %>%
- filter(Group %in% c(creeks, "Bee", "PID")) %>%
- filter(Odor %in% c("2-heptanone"), !Frequency %in% 10, !PulseDuration %in% 300)
- data <- droplevels(data)
- data$Odor <- factor(data$Odor, levels = levels(data$Odor),
- labels = c("2-heptanone"))
- # calculate the mean over all odor-pulseduration configurations
- d_mean <- data %>%
- group_by(Group, Winged, Odor, PulseDuration) %>%
- summarise_at(.funs = median, .vars=paste0("t", 1:5000))
- # calculate sample size (n) over all odor-pulseduration configurations
- d_n <- data %>%
- group_by(Group, Winged, Odor, PulseDuration) %>%
- summarise(N=n()) %>%
- group_by(Group, Winged) %>%
- summarise(N=unique(N))
- d_max <- d_mean %>%
- group_by(Group, Odor, PulseDuration) %>%
- summarise(Max=max(t400))
- # create d_all and add the sample size
- d_all <- d_mean
- d_all$PulseDuration <- factor(d_all$PulseDuration, levels = c('15', '30', '150'))
- d_all$PulseDur <- factor(d_all$PulseDuration)
- data_long <- gather(d_all, time, Voltage, t1:t5000, factor_key=TRUE)
- data_long$time <- (as.numeric(gsub("t", "", data_long$time)) - 200)/1000
- print(d_n)
- geom.text.size <- 5.2
- data_long$Winged <- factor(data_long$Winged,
- levels = levels(data_long$Winged),
- labels = c("Full-winged", "Wing-reduced", "none"))
- data_segm1 <- data.frame(y=c(0, 0), yend=c(0.1, 0), x=c(Inf, 2), xend=c(Inf, 3),
- PulseDuration = factor(150, levels = c(15, 30, 150)),
- Group=factor(c("Lug", 'Mt Burns'), levels = creeks),
- Odor=c("2-heptanone"))
- label1 <- data.frame(x=c(Inf, 3), y=c(0.05, 0.02),
- text=c("0.1 mV", '1 s'),
- PulseDuration= factor(150, levels = c(15, 30, 150)),
- Group=factor(c("Lug", 'Mt Burns'), levels = creeks),
- Odor=c("2-heptanone"))
- fw_text <- data.frame(x = 0, y = 0.02,
- text=c("full-winged"),
- PulseDuration= factor(15, levels = c(15, 30, 150)),
- Group=factor("Mt Burns", levels = creeks),
- Odor=c("2-heptanone"))
- wr_text <- data.frame(x = 0, y = 0.05,
- text=c("wing-reduced"),
- PulseDuration = factor(15, levels = c(15, 30, 150)),
- Group=factor("Mt Burns", levels = creeks),
- Odor=c("2-heptanone"))
- pulseDur_text <- data.frame(x = 0.25, y = 0.1,
- text=c("15 ms", "30 ms", "150 ms"),
- PulseDuration = factor(c("15", "30", "150"), levels = c("15", "30", "150")),
- Group=factor("Black Jacks", levels = creeks),
- Odor=c("2-heptanone"))
- data_long$Group <- factor(data_long$Group, levels = c(creeks, 'Bee', 'PID'))
- col <- c(fw_l, wr_l)
- Fen <- ggplot(data_long %>% filter(Group %in% creeks), aes(x=time, y=Voltage)) +
- facet_grid(PulseDuration ~ Group, switch = "y") +
- geom_rect(mapping=aes(xmin=0, xmax=as.numeric(as.character(PulseDuration))/1000,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_line(aes(colour = Winged)) +
- theme_bw() +
- labs(y = ' ') +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_blank(),
- axis.ticks = element_blank(),
- axis.title.x = element_blank(),
- panel.border = element_blank(),
- plot.margin = unit(c(5.5, 5.5, 5.5, 5.5), "points"),
- strip.background = element_rect(colour="NA", fill="NA"),
- text = element_text(size=18),
- aspect.ratio = ar * 3,
- strip.text.y.left = element_blank(),
- strip.text.x = element_text(size = 20),
- axis.text=element_blank()) +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- geom_text(data = wr_text, aes(x=x, y=y, label=text), size=6, hjust=0, color=wr_l) +
- geom_text(data = fw_text, aes(x=x, y=y, label=text), size=6, hjust=0, color=fw_l) +
- geom_text(data = pulseDur_text, aes(x=x, y=y, label=text), size=5, hjust=0, color='black') +
- # geom_text(data = Diffs, aes(x=0.4, y=Max + 0.025, label=Star), size=geom.text.size, hjust=0.5) +
- coord_cartesian(expand = TRUE,
- clip = 'off')
- Fen
- col <- c('blue')
- data_segm1 <- data.frame(y=0, yend=1, x=Inf, xend=Inf,
- PulseDuration = factor(150, levels = c(15, 30, 150)),
- Group=factor("Bee"),
- Odor=c("2-heptanone"))
- label1 <- data.frame(x=Inf, y=0.5,
- text=c("1 mV"),
- PulseDuration= factor(150, levels = c(15, 30, 150)),
- Group=factor("Bee"),
- Odor=c("2-heptanone"))
- Bee <- ggplot(data_long %>% filter(Group %in% 'Bee'), aes(x=time, y=Voltage)) +
- facet_grid(PulseDuration ~ Group, switch = "y") +
- geom_rect(mapping=aes(xmin=0, xmax=as.numeric(as.character(PulseDuration))/1000,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_line(aes(colour = Winged)) +
- theme_bw() +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_blank(),
- axis.title = element_blank(),
- plot.margin = unit(c(5.5, 5.5, 5.5, 5.5), "points"),
- strip.background = element_rect(colour="NA", fill="NA"),
- panel.border = element_blank(),
- text = element_text(size=18),
- axis.ticks = element_blank(),
- aspect.ratio = ar * 3,
- strip.text.y = element_blank(),
- strip.text.x = element_text(size = 20),
- axis.text=element_blank()) +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- # geom_text(data = Diffs, aes(x=0.4, y=Max + 0.025, label=Star), size=geom.text.size, hjust=0.5) +
- coord_cartesian(expand = TRUE,
- clip = 'off')
- Bee
- data_segm1 <- data.frame(y=0, yend=3, x=Inf, xend=Inf,
- PulseDuration = factor(150, levels = c(15, 30, 150)),
- Group=factor("PID"),
- Odor=c("2-heptanone"))
- label1 <- data.frame(x=Inf, y=1.5,
- text=c("3 V"),
- PulseDuration= factor(150, levels = c(15, 30, 150)),
- Group=factor("PID"),
- Odor=c("2-heptanone"))
- PID <- ggplot(data_long %>% filter(Group %in% 'PID'), aes(x=time, y=Voltage)) +
- facet_grid(PulseDuration ~ Group, switch = "y") +
- geom_rect(mapping=aes(xmin=0, xmax=as.numeric(as.character(PulseDuration))/1000,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_line(aes(colour = Winged)) +
- theme_bw() +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line = element_blank(),
- axis.title = element_blank(),
- plot.margin = unit(c(5.5, 5.5, 5.5, 5.5), "points"),
- strip.background = element_rect(colour="NA", fill="NA"),
- panel.border = element_blank(),
- text = element_text(size=18),
- axis.ticks = element_blank(),
- aspect.ratio = ar * 3,
- strip.text.y = element_blank(),
- strip.text.x = element_text(size = 20),
- axis.text=element_blank()) +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- # geom_text(data = Diffs, aes(x=0.4, y=Max + 0.025, label=Star), size=geom.text.size, hjust=0.5) +
- coord_cartesian(expand = TRUE,
- clip = 'off')
- PID
- FigS1A <- egg::ggarrange(Fen, Bee, PID,
- heights = 5, nrow = 1, ncol=3, widths = c(5, 1, 1), debug = FALSE)
- ########################
- #### Function to plot parameters for Response Dynamic:
- plot_bayes_time <- function(data, creeks, yvalue, ylab) {
-
- data <- data2
-
- data <- data %>%
- filter(PulseDuration==300, Group %in% creeks, Odor %in% odors)
-
- df <- as.data.frame(matrix(NA, ncol=7, nrow = 0))
- df_diff <- as.data.frame(matrix(NA, ncol=8, nrow = 0))
-
- ##########################################################################
-
- for (cc in seq_along(creeks)) {
-
- subdata <- data %>%
- filter(Group == creeks[cc])
-
- subdata$PulseDuration <- as.factor(subdata$PulseDuration)
- subdata <- droplevels(subdata)
-
- m1 <- lm(log2(get(yvalue)) ~ Winged * Odor, data = subdata)
- par(mfrow=c(2,2))
- plot(m1)
- summary(m1)
-
- subdata$Odor
-
- newdat <- expand.grid(Winged = factor(c("full-winged","wing-reduced")),
- PulseDuration = factor(levels(subdata$PulseDuration), levels = levels(subdata$PulseDuration)),
- Odor = factor(levels(subdata$Odor), levels = odors))
-
- bsim <- sim(m1, n.sim=nsim)
-
- fitmat <- matrix(ncol=nsim, nrow=nrow(newdat))
- Xmat <- model.matrix( ~ Winged * Odor, data=newdat)
- for(i in 1:nsim) fitmat[,i] <- Xmat %*% bsim@coef[i,] #fitted values
-
- tmp <- t(apply(fitmat, 1, quantile, prob=c(0.025, 0.5, 0.975)))
-
- tst <- cbind(newdat, tmp, 'Group' = creeks[cc])
-
- df <- rbind(df, tst)
-
- for (oo in odors) {
-
- dat <- fitmat[which(newdat$PulseDuration == 300 & newdat$Odor == oo), ]
- star2 <- sum(dat[1, ] > dat[2, ]) / nsim
-
-
-
- df_diff <- rbind(df_diff,
- data.frame('Group' = creeks[cc],
- 'Odor' = oo,
- 'Diff' = star2))
-
- }
-
-
- }
-
- df_diff
- df_diff$Star1 <- df_diff$Diff >= 0.95 | df_diff$Diff <= 0.05
- df_diff$Star2 <- df_diff$Diff >= 0.99 | df_diff$Diff <= 0.01
-
-
- colnames(df)[4:6] <- c('lower', 'fitted', 'upper')
-
-
- df$lower <- 2^(df$lower)
- df$fitted <- 2^(df$fitted)
- df$upper <- 2^(df$upper)
-
-
- subdf <- df %>%
- group_by(Group, Odor) %>%
- summarise(max=max(fitted))
-
-
- sub_data2 <- data %>%
- filter(Group!="Fenestrata-Bee",
- Group!="Bee",
- Odor!="Blank",
- Odor!="Propionic acid",
- Odor!="Stonefly")
-
- sub_data2 <- droplevels(sub_data2)
- df <- droplevels(df)
- df_diff <- droplevels(df_diff)
-
- sub_data2$PulseDuration <- factor(sub_data2$PulseDuration,
- levels = c(300))
-
- df$PulseDuration <- factor(df$PulseDuration,
- levels = c(300))
-
-
- df$NewGroup <- as.factor(paste(df$Winged, df$Group))
-
- df_diff$NewGroup <- as.factor(paste(df_diff$Winged, df_diff$Group))
-
-
-
- col <- c(fw_l, wr_l)
- col <- c(col[1], "black", col[2])
-
- df_diff$Star[df_diff$Star1 == TRUE] <- "*"
- df_diff$Star[df_diff$Star1 == FALSE] <- " "
- df_diff$Star[df_diff$Star2 == TRUE] <- "**"
- stars <- df_diff$Star
-
- df_diff$NewGroup <- "star"
-
- print(df_diff)
-
- df_diff$Winged <- "star"
-
- df$Odor <- factor(df$Odor, levels = odors)
- df_diff$Odor <- factor(df_diff$Odor, levels = odors)
- sub_data2$Odor <- factor(sub_data2$Odor, levels = odors)
-
- MAX <- sub_data2 %>%
- group_by(Odor, Group) %>%
- summarise(max = max(get(yvalue), na.rm = TRUE))
- df_diff <- merge(df_diff, MAX, by=c('Odor', 'Group'))
-
- p <- ggplot(df, aes(x=Group, y=fitted, color=Winged)) +
- geom_hline(aes(yintercept=0), color="black", size=0.3, linetype="dotted") +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- scale_y_continuous(breaks = pretty_breaks(n=2), expand = expansion(mult = c(0.05, .1))) +
- facet_grid(Odor ~ Group, scales = "free") +
- geom_beeswarm(sub_data2, dodge.width = 0.8, mapping=aes(x=Group, y=get(yvalue)), cex=1, shape=21) +
- geom_point(aes(group=NewGroup), position = position_dodge(0.8), shape="_", cex=4, color="black") +
- geom_errorbar(aes(ymin=lower, ymax=upper, group=Winged),
- color="black", width=0, position=position_dodge(0.8), size=0.5) +
- theme_bw() +
- ylab(ylab) +
- xlab(" ") +
- scale_x_discrete(labels=c("", "")) +
- geom_text(data=df_diff, aes(x=Group, y = max * 0.9, label=Star), size = 7) +
- theme(legend.position = "none",
- panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- axis.line.x = element_blank(),
- axis.line.y = element_line(colour = "black"),
- axis.title.x=element_blank(),
- axis.ticks.x=element_blank(),
- strip.text.x = element_text(size = 16),
- strip.background = element_rect(colour="NA", fill="NA"),
- strip.text.y = element_text(size = 16, angle = 0, hjust = 0),
- panel.border = element_blank(),
- text = element_text(size=18),
- axis.text=element_text(colour="black"))
- print(p)
-
- return(p)
- }
- ####### Figure 2 ##########################
- Fig2A <- plot_bayes_time(data = data2, yvalue = "Timeto10", ylab = "Response onset (ms)", creeks = creeks)
- Fig2B <- plot_bayes_time(data2, "TimeBackTo10", "Response offset (ms)", creeks = creeks)
- ################## 10 Hz Stimulations ###########
- on_times <- seq(0, 2.9, by = 0.1)
- off_times <- seq(0.05, 3, by=0.1)
- times_tst_on <- expand.grid(levels(data2$Odor), "on"=on_times)
- times_tst_off <- expand.grid(levels(data2$Odor), "off"=off_times)
- Times <- data.frame("Odor"=times_tst_on$Var1, "on"=times_tst_on$on, "off"=times_tst_off$off)
- Times <- Times[!(Times$Odor=="Blank"), ]
- plot_freq <- function(data, group, times) {
-
- data <- data %>%
- filter(PulseDuration == 50,
- Odor == "2-heptanone")
-
- data <- droplevels(data)
-
- d_mean <- data %>%
- filter(Group %in% group) %>%
- group_by(Group, Winged, Odor, PulseDuration) %>%
- summarise_at(.funs = median, .vars=paste0("t", 1:5000))
-
- # calculate the standard error
- d_se <- data %>%
- filter(Group %in% group) %>%
- group_by(Group, Winged, Odor, PulseDuration) %>%
- summarise_at(.funs = se, .vars=paste0("t", 1:5000))
-
- # calculate sample size (n)
- d_n <- data %>%
- filter(Group %in% group) %>%
- group_by(Group, Winged, Odor, PulseDuration) %>%
- summarise(N=n()) %>%
- group_by(Group, Winged) %>%
- summarise(N=unique(N))
- d_n <- unite(d_n, col = "Winged", sep = " (n = ")
- d_n$t <- ")"
- d_n <- unite(d_n, col = "Winged", sep = "")
-
-
- d_all <- d_mean
-
-
- data_se <- gather(d_se, time, se, t1:t5000, factor_key=TRUE)
-
- d_all$Odor <- factor(d_all$Odor)
- times <- times[times$Odor == "2-heptanone", ]
-
- data_long <- gather(d_all, time, Voltage, t1:t5000, factor_key=TRUE)
- data_long$time <- (as.numeric(gsub("t", "", data_long$time)) - 200)/1000
- data_long$SE <- data_se$se
- data_long$PulseDuration <- factor(data_long$PulseDuration, levels = c("50"), labels = c("10 Hz"))
-
-
- col <- c("grey", "#00DC00")
-
- data_segm1 <- data.frame(y=c(0, -0.02), yend=c(0.05, -0.02), x=c(Inf, 1), xend=c(Inf, 2),
- PulseDuration = factor(150, levels = c(15, 30, 150)),
- Group=factor(c("Six Mile", 'Mt Burns')),
- Odor=c("2-heptanone"))
- label1 <- data.frame(x=c(Inf, 2), y=c(0.035, -0.035),
- text=c("0.05 mV", '1 s'),
- PulseDuration= factor(150, levels = c(15, 30, 150)),
- Group=factor(c("Six Mile", 'Mt Burns')),
- Odor=c("2-heptanone"))
-
- ## For Lug, we stimultaed with a different stimulus protocol were the stimulation
- ## with several pulses in a specific frequency was not 10 Hz and the total length
- ## o f stimulation was less than 3 seconds. We exclude these data here.
- dat <- data_long %>% filter(Group == creeks, Group != 'Lug')
- dat$Group <- factor(dat$Group, levels = c(creeks))
-
- W <- ggplot(dat, aes(x=time, y=Voltage)) +
- geom_rect(data=times, inherit.aes=FALSE, fill = "grey95",
- aes(xmin=on, xmax=off, ymin=-Inf, ymax=Inf)) +
- # geom_ribbon(aes(ymin=Voltage + SE, ymax=Voltage - SE, fill=Winged), color=NA, alpha=0.3) +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_line(aes(colour = Winged)) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- facet_grid(. ~ Group, scales = "free_y", drop = FALSE) +
- scale_y_continuous(breaks = pretty_breaks(n = 3)) +
- theme_bw() +
- labs(y = " ") +
- guides(fill = 'none') +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- coord_cartesian(xlim =c(0, 3.5)) +
- theme(panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- legend.position = "none",
- strip.background = element_rect(colour="NA", fill="NA"),
- strip.text.y = element_blank(),
- legend.title = element_blank(),
- axis.line = element_blank(),
- axis.ticks = element_blank(),
- axis.title.x = element_blank(),
- axis.title.y = element_text(size = 35),
- panel.border = element_blank(),
- axis.text=element_blank(),
- strip.text = element_blank(),
- text = element_text(size=18))
- W
-
- col <- c("blue")
- data_segm1 <- data.frame(y=0, yend=1, x=Inf, xend=Inf,
- PulseDuration = factor(150, levels = c(15, 30, 150)),
- Group=factor("Bee"),
- Odor=c("2-heptanone"))
- label1 <- data.frame(x=Inf, y=0.5,
- text=c("1 mV"),
- PulseDuration= factor(150, levels = c(15, 30, 150)),
- Group=factor("Bee"),
- Odor=c("2-heptanone"))
-
- B <- ggplot(data_long[data_long$Group %in% 'Bee', ], aes(x=time, y=Voltage)) +
- geom_rect(data=times, inherit.aes=FALSE, fill = "grey95",
- aes(xmin=on, xmax=off, ymin=-Inf,ymax=Inf)) +
- # geom_ribbon(aes(ymin=Voltage + SE, ymax=Voltage - SE, fill=Winged), color=NA, alpha=0.3) +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_line(aes(colour = Winged)) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- facet_grid(Odor ~ Group, scales = "free_y") +
- scale_y_continuous(breaks = pretty_breaks(n = 3)) +
- theme_bw() +
- coord_cartesian(xlim =c(0, 3.5)) +
- # guides(fill = 'none') +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- theme(panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- legend.position = "none",
- strip.background = element_rect(colour="NA", fill="NA"),
- panel.border = element_blank(),
- legend.title = element_blank(),
- axis.line = element_blank(),
- strip.text = element_blank(),
- axis.text=element_blank(),
- axis.ticks = element_blank(),
- axis.title = element_blank(),
- text = element_blank())
- B
-
-
- data_segm1 <- data.frame(y=0, yend=1, x=3, xend=3,
- PulseDuration = factor(150, levels = c(15, 30, 150)),
- Group=factor("PID"),
- Odor=c("2-heptanone"))
- label1 <- data.frame(x=3, y=0.5,
- text=c("1 V"),
- PulseDuration= factor(150, levels = c(15, 30, 150)),
- Group=factor("PID"),
- Odor=c("2-heptanone"))
-
- P <- ggplot(data_long[data_long$Group %in% 'PID', ], aes(x=time, y=Voltage)) +
- geom_rect(data = times, inherit.aes = FALSE, fill = "grey95",
- aes(xmin = on, xmax=off, ymin = -Inf, ymax = Inf)) +
- # geom_ribbon(aes(ymin=Voltage + SE, ymax=Voltage - SE, fill=Winged), color=NA, alpha=0.3) +
- geom_hline(mapping = aes(yintercept = 0), color="black", size=0.3, linetype="dotted") +
- geom_line(aes(colour = Winged)) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- facet_grid(Odor ~ Group, scales = "free_y") +
- scale_y_continuous(breaks = pretty_breaks(n = 3)) +
- theme_bw() +
- coord_cartesian(xlim =c(0, 3.5)) +
- geom_segment(data=data_segm1,
- aes(x=x, y=y, yend=yend, xend=xend)) +
- geom_text(data = label1, aes(x=x, y=y, label=text), size=geom.text.size, hjust=1.1) +
- # guides(fill = 'none') +
- theme(panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- legend.position = "none",
- strip.background = element_rect(colour="NA", fill="NA"),
- panel.border = element_blank(),
- legend.title = element_blank(),
- axis.ticks = element_blank(),
- axis.text=element_blank(),
- axis.line = element_blank(),
- axis.title = element_blank(),
- strip.text = element_blank(),
- text = element_blank())
- P
-
- list(W = W, B = B, P = P)
-
- }
- FigS1C_new <- plot_freq(data2, c(creeks, 'Bee', 'PID'), Times)
- ## Power spectrum density analysis ####
- ## whole frequency band mean over groups of individuals ##
- ## Freq = 10 for individual ID ##
- spec <- data.frame(x=as.numeric(NULL), y=as.numeric(NULL), sd=as.numeric(NULL), se=as.numeric(NULL),
- Odor=as.character(NULL), Winged=as.character(NULL),
- Group=as.character(NULL), TimeWindow=as.character(NULL))
- spec_ind <- data.frame(y=as.numeric(NULL),
- yMean=as.numeric(NULL),
- xMean=as.numeric(NULL),
- zMean=as.numeric(NULL),
- yMax=as.numeric(NULL),
- xMax=as.numeric(NULL),
- zMax=as.numeric(NULL),
- id=as.numeric(NULL),
- Odor=as.character(NULL), Winged=as.character(NULL),
- Group=as.character(NULL), TimeWindow=as.character(NULL))
- winged <- as.vector(t(unique(data2[, "Winged"])))
- for(k in 1:length(winged)){
- group <- as.vector(t(unique(data2[data2$Winged==winged[k], "Group"])))
-
- for (l in 1:length(group)) {
- odors <- as.vector(t(unique(data2[data2$Winged==winged[k] & data2$Group==group[l], "Odor"])))
-
- for (i in 1:length(odors)){
- Hept <- data2 %>%
- filter(Odor==odors[i], Winged==winged[k], Group==group[l], (AntennaState=="alive" | is.na(AntennaState)))
-
-
- ### full length of 3 seconds ####
- spect <- apply(Hept[,paste0("t", 601:3600)],
- MARGIN = 1,
- FUN=function(x) spec.mtm(ts(x, frequency = 1000), taper = "sine", plot=F))
-
- spec <- rbind(spec, data.frame(x=spect[[1]]$freq[1:329],
- y=rowMeans(sapply(spect, '[[', 4))[1:329],
- sd=apply((sapply(spect, '[[', 4)),1, sd)[1:329],
- se=apply((sapply(spect, '[[', 4)),1, se)[1:329],
- Odor=odors[i],
- Winged=winged[k],
- Group=group[l],
- TimeWindow="T0-3"))
-
- # select frequencies between 9.5 and 10.5 and calculate the mean values
- sel_row_10 <- which(spect[[1]]$freq > 9 & spect[[1]]$freq < 11)
- sel_row_pre10 <- which(spect[[1]]$freq > 8 & spect[[1]]$freq < 9)
- sel_row_post10 <- which(spect[[1]]$freq > 11 & spect[[1]]$freq < 12)
-
- spec_ind <- rbind(spec_ind, data.frame(y=colSums(sapply(spect, '[[', 4)[sel_row_10,]),
- yMean=colMeans(sapply(spect, '[[', 4)[sel_row_10,]),
- xMean=colMeans(sapply(spect, '[[', 4)[sel_row_pre10,]),
- zMean=colMeans(sapply(spect, '[[', 4)[sel_row_post10,]),
- yMax=apply(sapply(spect, '[[', 4)[sel_row_10,], 2, max),
- xMax=apply(sapply(spect, '[[', 4)[sel_row_pre10,], 2, min),
- zMax=apply(sapply(spect, '[[', 4)[sel_row_post10,], 2, min),
- id=as.character(Hept$ID_new),
- Odor=odors[i],
- Winged=winged[k],
- Group=group[l],
- TimeWindow="T0-3"))
-
- }
- }
-
- }
- spec_sel <- spec %>%
- filter(Odor!="Blank", Odor!="Stonefly", Odor!="Propionic acid", TimeWindow=="T0-3")
- spec_sel <- droplevels(spec_sel)
- spec_sel$Odor <- factor(spec_sel$Odor, levels = c("2-heptanone",
- "1-octanol",
- "2-butanone"))
- col <- c("#FF00FF", "#00DC00")
- plot_freq <- function(data) {
-
- data <- data %>%
- filter(Odor=="2-heptanone")
-
- col <- c(fw_l, wr_l, "black")
-
- data$TimeWindow <- factor(data$TimeWindow, levels = levels(data$TimeWindow),
- labels = "")
-
- odors <- data.frame("Odor"=levels(data$Odor))
-
- ## For Lug, we stimultaed with a different stimulus protocol were the stimulation
- ## with several pulses in a specific frequency was not 10 Hz and the total length
- ## o f stimulation was less than 3 seconds. We exclude these data here.
- dat <- data %>% filter(Group %in% creeks & Group != 'Lug')
- dat$Group <- factor(dat$Group, levels = c(creeks))
-
- W <- ggplot(dat,
- aes(x=x, y=log10(y))) +
- geom_rect(mapping=aes(xmin=9.5, xmax=10.5,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- geom_ribbon(aes(ymin=log10(y + se), ymax=log10(y - se), fill=Winged),
- color=NA, alpha=0.3) +
- geom_line(aes(colour = Winged)) +
- facet_grid(. ~ Group, drop = FALSE) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = col) +
- theme_bw() +
- scale_y_continuous(breaks = pretty_breaks(n=3)) +
- theme(panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- legend.position = "none",
- text = element_text(size=18),
- axis.title = element_blank(),
- axis.line = element_line(colour = "black"),
- strip.background = element_rect(colour="NA", fill="NA"),
- strip.text = element_blank(),
- panel.border = element_blank(),
- axis.text=element_text(colour="black"))
- W
-
- B <- ggplot(data %>% filter(Group %in% c("Bee")),
- aes(x=x, y=log10(y))) +
- geom_rect(mapping=aes(xmin=9.5, xmax=10.5,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- geom_ribbon(aes(ymin=log10(y + se), ymax=log10(y - se), fill=Winged),
- color=NA, alpha=0.3) +
- geom_line(aes(colour = Winged)) +
- facet_grid(. ~ Group) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = "blue") +
- theme_bw() +
- scale_y_continuous(breaks = pretty_breaks(n=3)) +
- theme(panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- legend.position = "none",
- text = element_text(size=18),
- axis.title = element_blank(),
- axis.line = element_line(colour = "black"),
- strip.background = element_rect(colour="NA", fill="NA"),
- strip.text = element_blank(),
- panel.border = element_blank(),
- axis.text=element_text(colour="black"))
- B
-
- P <- ggplot(data %>% filter(Group %in% c("PID")),
- aes(x=x, y=log10(y))) +
- geom_rect(mapping=aes(xmin=9.5, xmax=10.5,
- ymin=-Inf, ymax=Inf),
- fill="grey90", alpha=1) +
- geom_ribbon(aes(ymin=log10(y + se), ymax=log10(y - se), fill=Winged),
- color=NA, alpha=0.3) +
- geom_line(aes(colour = Winged)) +
- facet_grid(. ~ Group) +
- scale_color_manual(aesthetics = c("colour", "fill"), values = "blue") +
- theme_bw() +
- scale_y_continuous(breaks = pretty_breaks(n=3)) +
- theme(panel.grid.major = element_blank(),
- panel.grid.minor = element_blank(),
- legend.position = "none",
- axis.title = element_blank(),
- text = element_text(size=18),
- axis.line = element_line(colour = "black"),
- strip.background = element_rect(colour="NA", fill="NA"),
- strip.text = element_blank(),
- panel.border = element_blank(),
- axis.text=element_text(colour="black"))
- P
-
- list1 <- list("traceW"=W,
- "traceB"=B,
- "traceP"=P)
- return(list1)
-
- }
- p5 <- plot_freq(data = spec_sel); p5
- #create common x and y labels
- FigS1D <- egg::ggarrange(p5$traceW, p5$traceB, p5$traceP,
- bottom = textGrob("Frequency (Hz)", gp = gpar(fontsize=16), vjust = 0.5),
- left = textGrob(expression(paste("Power (", mV^2/Hz, ") [log10]")),
- gp = gpar(fontsize=18), vjust = 0.5, rot = 90),
- ncol = 3, nrow=1,
- label.args = list(gp=gpar(fontface="bold", fontsize=18), hjust=0, vjust=1),
- widths = c(5,1,1),
- debug = FALSE)
- FigS1D
- FigS1C <- egg::ggarrange(FigS1C_new$W, FigS1C_new$B, FigS1C_new$P,
- ncol = 3, nrow=1, heights = 1,
- widths = c(5, 1, 1),
- debug = FALSE)
- FigS1C
- Fig1AB <- cowplot::plot_grid(Fig1A, Fig1B,
- align = "hv", axis = "trbl",
- labels = c("A", "B"),
- rel_widths = c(1.2, 1.5),
- nrow = 1, label_size = 18)
- ### Fig. 1
- cowplot::plot_grid(Fig1AB, Fig1C, nrow = 2, rel_heights = c(1, 1.3),
- labels = c("", "C"), label_size = 18)
- ggsave("Figure1.pdf", path = fig_path, width = 12, height = 10)
- ### Fig. 2
- cowplot::plot_grid(Fig2A, Fig2B, nrow = 2,
- labels = LETTERS[1:2], label_size = 18)
- ggsave("Figure2.pdf", path = fig_path, width = 9, height = 12)
- ### Fig. S1
- cowplot::plot_grid(FigS1A, FigS1B, NULL, FigS1C, FigS1D, nrow = 5, rel_heights = c(1.2, 1, 0.05, 0.3, 0.8),
- labels = c('A', 'B', '', 'C', 'D'), label_size = 18)
- ggsave("FigureS1.pdf", path = fig_path, width = 14, height = 14)
|