Denecke et. al 2017: Describing the Role of Drosophila melanogaster ABCB Transporters in Insecticide Biology Using CRISPR-Cas9 Knockouts
Prior To Beginning
Welcome to the R markdown walk through for the paper Welcome to the R markdown companion to the recent publication Describing the Role of Drosophila melanogaster ABCB Transporters in Insecticide Biology Using CRISPR-Cas9 Knockouts (Denecke et. al 2017). This document is designed to guide the reader through every step of the analysis accomplished in the paper and allow for the reproducibility of the study.
Before beginning a few housekeeping issues must be addressed. Firstly, all code demonstrated here was optimized on a Linux machine but has been tested on a mac and a PC. While it should work the particular specifications of the environment in which it was written are listed below.
version## _
## platform x86_64-redhat-linux-gnu
## arch x86_64
## os linux-gnu
## system x86_64, linux-gnu
## status
## major 3
## minor 4.1
## year 2017
## month 06
## day 30
## svn rev 72865
## language R
## version.string R version 3.4.1 (2017-06-30)
## nickname Single CandleNext, the analysis will require several packages. These packages can be installed and loaded into the R environment with the following chunk of code.
list.of.packages <- c("devtools","plotrix")
new.packages <- list.of.packages[!(list.of.packages %in% installed.packages()[,"Package"])]
if(length(new.packages)) install.packages(new.packages)
if(!"insect.toxicology" %in% installed.packages()) devtools::install_github("shanedenecke/insect.toxicology")
library(insect.toxicology)
library(plotrix)The data will need to be downloaded from the Mendeley repository. This can be done with the following code.
download.file("https://github.com/shanedenecke/Denecke_2017_ABCB_Insecticide/raw/master/ABCB_raw.data.binary",destfile="ABCB_raw.data.binary")
load("~/test_ABCB/ABCB_raw.data.binary")Lastly, some custom functions that are used during the analysis can be loaded now.
firstup <- function(x) {
substr(x, 1, 1) <- toupper(substr(x, 1, 1))
x
}
multi.gsub <- function(vector,find,replace){
for (i in 1:length(find)){
vector <- gsub(find[i],replace[i],vector)
}
return(vector)
} Figure 1
Now the analysis proper can be started. The first figure in the paper deals with the change in insecticide resistance engendered by each of the 3 ABCB deletions.
First import the toxicology data into R.
knitr::kable(head(Fig1_ABCB_Combined))| X | pesticide | genotype | rep | dose | alive | dead | total | percent_survival | percent_mortality | stage | date | conditions |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 1 | nitenpyram | wxac | 1 | 0.2 | 45 | 5 | 50 | 0.90 | 0.10 | L1 | 7.3.16 | 25d_Hoff |
| 2 | nitenpyram | wxac | 2 | 0.2 | 41 | 9 | 50 | 0.82 | 0.18 | L1 | 7.3.16 | 25d_Hoff |
| 3 | nitenpyram | wxac | 3 | 0.2 | 40 | 10 | 50 | 0.80 | 0.20 | L1 | 7.3.16 | 25d_Hoff |
| 4 | nitenpyram | wxac | 4 | 0.2 | 45 | 5 | 50 | 0.90 | 0.10 | L1 | 7.3.16 | 25d_Hoff |
| 5 | nitenpyram | wxac | 5 | 0.2 | 42 | 8 | 50 | 0.84 | 0.16 | L1 | 7.3.16 | 25d_Hoff |
| 6 | nitenpyram | mdr49ko | 1 | 0.2 | 29 | 21 | 50 | 0.58 | 0.42 | L1 | 7.3.16 | 25d_Hoff |
The next chunk of code will add an abbots correction to for each insecticide/dose/genotype combination. This function (abbots.correct.table) is part of the insect.toxicology package available in my github repository (shanedenecke/insect.toxicology).
abcb.abbots <- abbots.correction.table(Fig1_ABCB_Combined,control="wxac",write=F) %>%
filter(dose!=0)This data must be formatted before plotting. We need to make calls about whether a genotype is more resistant or susceptible to a given genotype. As discussed in the paper, this is done by the rather crude method of non overlapping 95% confidence intervals outlined by Rosenheim and Hoy (1987). Once calls are made, they can be converted to numeric values for plotting.
index <- abcb.abbots %>% select(genotype,pesticide) %>% unique()
summary.matrix <- matrix(data=NA,nrow=length(unique(index$pesticide)),ncol=length(unique(index$genotype)),dimnames=list(unique(index$pesticide),unique(index$genotype)))
control="wxac"
for(i in 1:nrow(index)){
test.genotype <- as.character(index[i,1])
test.pesticide <- as.character(index[i,2])
sub <- subset(abcb.abbots, (genotype==test.genotype | genotype==control) & pesticide==test.pesticide)
emvec <- c()
for(d in unique(sub$dose)){
sub2 <- subset(sub,dose==d)
if(sub2$lcl[which(sub2$genotype=="wxac")] > sub2$ucl[which(sub2$genotype==test.genotype)]){
emvec <- c(emvec,"sus")
}else if (sub2$ucl[which(sub2$genotype=="wxac")] < sub2$lcl[which(sub2$genotype==test.genotype)]){
emvec <- c(emvec,"res")
}else{
emvec <- c(emvec,"nsd")
}
}
names(emvec) <- unique(sub$dose)
if(("res" %in% emvec) & ("sus" %in% emvec)){
summary.matrix[test.pesticide,test.genotype] <- 0
} else if("res" %in% emvec){
summary.matrix[test.pesticide,test.genotype] <- 1
} else if("sus" %in% emvec){
summary.matrix[test.pesticide,test.genotype] <- -1
} else {
summary.matrix[test.pesticide,test.genotype] <- 0
}
}
abcb.hm <- data.frame(summary.matrix,stringsAsFactors = F)
colnames(abcb.hm) <- multi.gsub(colnames(abcb.hm),find=c("wxac","mdr49ko","mdr50ko","mdr65ko"),c("Wxac","Mdr49KO","Mdr50KO","Mdr65KO"))
rownames(abcb.hm) <- sapply(rownames(abcb.hm),firstup)
abcb.hm$Wxac <- NULLNow, the data can be plotted using the plotrix package.
par(mar = c(5, 12, 4, 2),family="serif")
color2D.matplot(abcb.hm,
show.values = FALSE,
axes = FALSE,
xlab = "",
ylab = "",
vcex = 2,
vcol = "black",
extremes = c("red4", "green4"),
show.legend = F,
nslices=3)
axis(3, at = seq_len(ncol(abcb.hm)) - 0.5,
labels = names(abcb.hm), tick = FALSE, cex.axis = 1.75, face="bold")
axis(2, at = seq_len(nrow(abcb.hm)) -0.5,
labels = rev(rownames(abcb.hm)), tick = FALSE, las = 1, cex.axis = 1.5,face="bold")
Supplementary Table 3
The raw data used to make this figure can be found in Supplementary Table 3
knitr::kable(head(abcb.abbots))| genotype | dose | pesticide | corrected.survival | conf.int | ucl | lcl |
|---|---|---|---|---|---|---|
| wxac | 0.2 | nitenpyram | 0.9434797 | -0.0731556 | 1.0166353 | 0.8703241 |
| mdr49ko | 0.2 | nitenpyram | 0.9886993 | -0.1808891 | 1.1695883 | 0.8078102 |
| mdr50ko | 0.2 | nitenpyram | 0.9324040 | -0.0844708 | 1.0168748 | 0.8479331 |
| mdr65ko | 0.2 | nitenpyram | 0.8298523 | -0.0847912 | 0.9146435 | 0.7450611 |
| wxac | 0.6 | nitenpyram | 0.8947554 | -0.0812738 | 0.9760293 | 0.8134816 |
| mdr49ko | 0.6 | nitenpyram | 0.9134176 | -0.1795396 | 1.0929572 | 0.7338780 |
Supplementary Figure 3
In order to observe the effect of these knockouts on fitness, survival in the absence of insecticide can be compared between the knockouts and the control.
control <- subset(Fig1_ABCB_Combined,
(pesticide=="pyriprole" & dose==0) |
(pesticide=="spinosad" & dose==0) |
(pesticide=="ivermectin" & dose==0))
control$genotype <- factor(control$genotype,levels = c("wxac","mdr49ko","mdr50ko","mdr65ko"))
gp = ggplot(data =control, aes(x = genotype, y = alive))
gp = gp + geom_dotplot(aes(fill = genotype), binaxis = "y",
stackdir = "center", binwidth = 0.1, colour = "black",
position = position_dodge(width = 0.5), dotsize = 10)
gp = gp + stat_summary(geom = "errorbar", colour = "black", size = 0.6,
fun.data = mean_cl_normal, position = position_dodge(width = 0.5), width = 0.5)
gp = gp + ylab("Number Emerging\n")
gp = gp + xlab("\nGenotype")
gp = gp + stat_summary(fun.y = "mean", fun.ymax = "mean",
fun.ymin = "mean", geom = "crossbar", colour = "black",
position = position_dodge(width = 0.5), width = 0.5)
gp = gp + scale_y_continuous(limits = c(0, 55))
gp <- gp + scale_fill_manual(values = c("grey50","red","blue","green"))
gp = gp + theme(text = element_text(size = 18, face = "bold"),
axis.text.x = element_text(angle = 60, hjust = 1,
size = 14), axis.text.y = element_text(size = 14,
family = "", face = "bold"), panel.background = element_rect(fill = NA),
panel.border = element_rect(colour = "black", fill = NA),
panel.grid.minor = element_blank(), panel.grid.major.x = element_blank(),
plot.title = element_text(hjust = 0.5))
print(gp)
It appears that only Mdr49KO shows increased mortality in the absence of insecticides.
mod <- with(control,aov(alive~genotype))
TukeyHSD(x=mod, which='genotype', conf.level=0.95)$genotype[,"p adj"]## mdr49ko-wxac mdr50ko-wxac mdr65ko-wxac mdr50ko-mdr49ko
## 0.001700731 0.946214373 0.990000004 0.008917079
## mdr65ko-mdr49ko mdr65ko-mdr50ko
## 0.004355571 0.994812131Supplementary Figure 4
Supplementary figures 1 and 2 are not data driven and are not discussed in this file. Supplementary figure 4, however, reflects the increased susceptibility to imidacloprid metabolites observed in Mdr65KO.
First, the data is imported.
head(FigS4_M65_Metabolites)## X pesticide genotype rep dose alive dead total percent_survival
## 1 156 imi-olefin mdr65ko 1 0.92 1 49 50 0.02
## 2 157 imi-olefin mdr65ko 2 0.92 1 49 50 0.02
## 3 158 imi-olefin mdr65ko 3 0.92 0 50 50 0.00
## 4 159 imi-olefin mdr65ko 4 0.92 0 50 50 0.00
## 5 160 imi-olefin mdr65ko 5 0.92 0 50 50 0.00
## 6 376 imi-olefin mdr65ko 1 0.80 31 19 50 0.62
## percent_mortality stage date conditions
## 1 0.98 1st_instar 05/05/17 25d_hoff
## 2 0.98 1st_instar 05/05/17 25d_hoff
## 3 1.00 1st_instar 05/05/17 25d_hoff
## 4 1.00 1st_instar 05/05/17 25d_hoff
## 5 1.00 1st_instar 05/05/17 25d_hoff
## 6 0.38 L1 7.3.16 25d_HoffThe data can be formatted so that the plot will look nice.
FigS4_M65_Metabolites.corrected <- abbots.correction.table(FigS4_M65_Metabolites,control="wxac",write=F)
annot <- data.frame(x = rep(-Inf, length(unique(FigS4_M65_Metabolites.corrected$genotype))), y = rep(Inf, length(unique(FigS4_M65_Metabolites.corrected$genotype))), pesticide=as.character(unique(FigS4_M65_Metabolites.corrected$pesticide)),labs=LETTERS[1:length(unique(FigS4_M65_Metabolites.corrected$genotype))])
FigS4_M65_Metabolites.corrected$genotype <- factor(FigS4_M65_Metabolites.corrected$genotype,levels=c("wxac","mdr65ko"))
FigS4_M65_Metabolites.corrected$pesticide <- factor(FigS4_M65_Metabolites.corrected$pesticide,levels=c("imi-olefin","imi-5-oh"))Lastly the Plot for supplementary figure 4. It appears like Mdr65KO is more susceptible to IMI-Olefin and IMI-5OH.
gp <- ggplot(FigS4_M65_Metabolites.corrected,aes(x=dose,y=corrected.survival,fill=genotype))
gp <- gp+facet_grid(.~pesticide,scales="free_x",labeller = label_parsed)
gp <- gp+geom_bar(position=position_dodge(width=.5),stat="identity",colour="black",width=.5)
gp <- gp+geom_errorbar(aes(ymin=corrected.survival+conf.int,ymax=corrected.survival-conf.int),
position=position_dodge(width=.5),stat="identity",width=.5)
gp <- gp+ylab("Corrected Survival \n")
gp <- gp+xlab("\n Concentration of Nitenpyram (ppm)")
gp <- gp+scale_fill_manual(values=c("grey25","skyblue"))
gp <- gp+scale_y_continuous(breaks=c(seq(0,1,.2)))
gp <- gp+geom_text(size=7,data=annot,aes(x, y, label=labs,fontface="bold",family="serif"), inherit.aes=FALSE,vjust=1.2,hjust=-.2)
gp <- gp+theme_bw(base_size=14,base_family="serif")
gp <- gp+theme(text=element_text(face="bold",family="serif"),panel.grid=element_blank(),axis.ticks.x=element_line(),
panel.border=element_rect(colour="black",fill=NA),strip.text=element_text(size=18,face="bold"),
strip.background=element_rect("grey95"),
axis.title=element_text(size=18),axis.text.x=element_text(angle=30,hjust=1))
print(gp) 
Figure 2
Figure 2 considers the Lc50 ratios of Mdr65KO to its background control.
The first chunk of code will import the data.
head(Fig2_RR_Summary)## Insecticide Lethal.Dose.Ratio LCL UCL
## 1 Chlorantraniliprole 1.123 0.9690 1.300
## 2 DDT 1.867 1.6980 2.052
## 3 Imidacloprid 1.184 1.0160 1.379
## 4 Ivermectin 1.385 1.1239 1.548
## 5 Nitenpyram 4.660 4.2370 5.138
## 6 Pyripole 7.956 7.1130 8.899Fig2_RR_Summary <- arrange(Fig2_RR_Summary,Lethal.Dose.Ratio)
Fig2_RR_Summary$Insecticide <- factor(Fig2_RR_Summary$Insecticide,levels=Fig2_RR_Summary$Insecticide)The data can then be plotted to reflect the estimated differences MDR65 substrate affinity.
gp <- ggplot(data=Fig2_RR_Summary,aes(x=Insecticide,y=Lethal.Dose.Ratio))
gp <- gp+geom_bar(position="dodge",stat="identity",colour="black",fill="grey50")
gp <- gp+geom_hline(yintercept = 1,linetype=2)
gp <- gp+geom_errorbar(aes(ymin=LCL,ymax=UCL),width=.7,size=.5)
gp <- gp+labs(y="Susceptibility Ratio\n (Mdr65KO/Wxac)\n",x="\nPesticide")
gp <- gp+scale_y_continuous(expand = c(0, 0),limits=c(0,9),
breaks=c(seq(0,8,2)),labels=c(seq(0,8,2)))
gp <- gp+geom_text(size=4,label="N.S",aes(1, 1.75,family="serif"))
gp <- gp+theme_bw(base_size=14,base_family="serif")
gp <- gp+theme(text=element_text(face="bold",family="serif"),panel.grid=element_blank(),axis.ticks.x=element_line(),
panel.border=element_rect(colour="black",fill=NA),
axis.title=element_text(size=18),axis.text.x=element_text(angle=30,hjust=1),legend.position = "none")
print(gp)
Figure 3
Verapamil is a ABCB inhibitor that acts as an insecticide synergist. Figure 3 assessed the impact of verapamil synergism with and without a functional copy of Mdr65.
The data for figure 3 is imported.
head(Fig3_Verpamil_Synergism)## Pesticide Synergist Synergist_Concentration Genotype Rep Dose Alive
## 1 Nitenpyram None 0 wxac 1 0.0 44
## 2 Nitenpyram None 0 wxac 2 0.0 48
## 3 Nitenpyram None 0 wxac 3 0.0 45
## 4 Nitenpyram None 0 wxac 4 0.0 44
## 5 Nitenpyram None 0 wxac 5 0.0 49
## 6 Nitenpyram None 0 wxac 1 1.2 22
## Dead Total Percent_Survival Percent_Mortality Stage Date Conditions
## 1 6 50 0.88 0.12 L1 27.4.17 25d_hoff
## 2 2 50 0.96 0.04 L1 27.4.17 25d_hoff
## 3 5 50 0.90 0.10 L1 27.4.17 25d_hoff
## 4 6 50 0.88 0.12 L1 27.4.17 25d_hoff
## 5 1 50 0.98 0.02 L1 27.4.17 25d_hoff
## 6 28 50 0.86 0.56 L1 27.4.17 25d_hoffThe data for figure 3 is formatted.
ver.corrected <- lapply(split(Fig3_Verpamil_Synergism,Fig3_Verpamil_Synergism$Synergist),abbots.correction.table,control="wxac",write=F)
ver.corrected[["None"]]$synergist <- "None"
ver.corrected[["Verapamil"]]$synergist <- "Verapamil"
ver.summary <- rbindlist(ver.corrected) %>% data.frame()
ver.summary$genotype <- factor(ver.summary$genotype,levels=c("wxac","Mdr65KO"))
for(i in 1:nrow(ver.summary)){
if(is.na(ver.summary[i,"genotype"])){ ver.summary[i,"genotype"] <- "Mdr65KO"}
}
annot <- data.frame(x = rep(-Inf, length(unique(ver.summary$genotype))), y = rep(Inf, length(unique(ver.summary$genotype))), genotype=as.character(unique(ver.summary$genotype)),labs=LETTERS[1:length(unique(ver.summary$genotype))])Figure 3 is plotted. It looks as though the only difference between verapamil + and verapamil - treatments is in the Wxac background, when Mdr65 is present.
gp <- ggplot(ver.summary,aes(x=dose,y=corrected.survival,fill=synergist,group=interaction(dose,synergist)))
gp <- gp+facet_grid(.~genotype,scales="free_x",labeller = label_parsed)
gp <- gp+geom_bar(position=position_dodge(width=.5),stat="identity",colour="black",width=.5)
gp <- gp+geom_errorbar(aes(ymin=corrected.survival+conf.int,ymax=corrected.survival-conf.int),
position=position_dodge(width=.5),stat="identity",width=.5)
gp <- gp+ylab("Corrected Survival \n")
gp <- gp+xlab("\n Concentration of Nitenpyram (ppm)")
gp <- gp+scale_fill_manual(values=c("grey50","red"))
gp <- gp+scale_y_continuous(breaks=c(seq(0,1,.2)))
gp <- gp+geom_text(size=7,data=annot,aes(x, y, label=labs,fontface="bold",family="serif"), inherit.aes=FALSE,vjust=1.2,hjust=-.2)
gp <- gp+theme_bw(base_size=14,base_family="serif")
gp <- gp+theme(text=element_text(face="bold",family="serif"),panel.grid=element_blank(),axis.ticks.x=element_line(),
panel.border=element_rect(colour="black",fill=NA),strip.text=element_text(size=18,face="bold"),
strip.background=element_rect("grey95"),
axis.title=element_text(size=18),axis.text.x=element_text(angle=30,hjust=1))
print(gp) 
Figure 4
Another way to confirm the role of ABCB transporters (in this case Mdr65) in insecticide resistance is to use genetic complementation. Here, the Mdr65KO genotype is crossed to genotypes which lack a region of their genomes that spans the Mdr65 locus.
The data for figure 4 is imported into R. Very minor formatting is also applied.
head(Fig4_Mdr65_Deficiency)## Pesticide Genotype Rep Dose Alive Dead Total Percent_Survival
## 1 Nitenpyram Mdr65KO x w1118 1 2 19 1 20 0.95
## 2 Nitenpyram Mdr65KO x w1118 2 2 20 0 20 1.00
## 3 Nitenpyram Mdr65KO x w1118 3 2 18 2 20 0.90
## 4 Nitenpyram Mdr65KO x w1118 4 2 20 0 20 1.00
## 5 Nitenpyram Mdr65KO x w1118 5 2 19 1 20 0.95
## 6 Nitenpyram Wxac x w1118 1 2 20 0 20 1.00
## Percent_Mortality Stage Date Conditions
## 1 0.05 Adult 23.9.16 22d
## 2 0.00 Adult 23.9.16 22d
## 3 0.10 Adult 23.9.16 22d
## 4 0.00 Adult 23.9.16 22d
## 5 0.05 Adult 23.9.16 22d
## 6 0.00 Adult 23.9.16 22dFig4_Mdr65_Deficiency$Genotype <- factor(Fig4_Mdr65_Deficiency$Genotype,levels=c("Wxac x w1118","Mdr65KO x w1118", "Mdr65KO x Excel_DF", "Mdr65KO x Bloom_DF", "Mdr65KO x Mdr65KO"))Figure 4 can then be plotted. The black and grey dots represent genotypes with 2 and 1 copy of Mdr65 respectively. These appear far more resistant than do the latter 3 genotypes which all lack Mdr65.
gp <- ggplot(Fig4_Mdr65_Deficiency,aes(x=Genotype,y=Alive,fill=Genotype))
gp <- gp+stat_summary(geom="errorbar",colour='black',size=1,width=.4,fun.data="mean_cl_normal")
gp <- gp+stat_summary(fun.y="mean",fun.ymax="mean",fun.ymin="mean", geom="crossbar", width=.4,colour='black',size=.8)
gp <- gp+geom_dotplot(binaxis='y',stackdir="center",binwidth=.75)
gp <- gp+scale_fill_manual(values=c("black","grey50","white","white","white"))
gp <- gp+ylab("Percent Survival\n")
gp <- gp+xlab("\nGenotype")
gp <- gp+theme_bw(base_size=14,base_family="serif")
gp <- gp+theme(text=element_text(face="bold",family="serif"),panel.grid=element_blank(),axis.ticks.x=element_line(),
panel.border=element_rect(colour="black",fill=NA),strip.text=element_text(size=18,face="bold"),
strip.background=element_rect("grey95"),
axis.title=element_text(size=18),axis.text.x=element_text(angle=30,hjust=1),legend.position = "none")
print(gp)
While these differences are pretty obvious, analysis of variance (ANOVA) with Tukey’s post-hoc honestly significant difference can be used to report p values for each pairwise comparison.
mod <- with(Fig4_Mdr65_Deficiency,aov(Alive~Genotype))
TukeyHSD(x=mod, which='Genotype', conf.level=0.95)$Genotype[,"p adj"]## Mdr65KO x w1118-Wxac x w1118
## 9.988116e-01
## Mdr65KO x Excel_DF-Wxac x w1118
## 4.973799e-14
## Mdr65KO x Bloom_DF-Wxac x w1118
## 2.275957e-14
## Mdr65KO x Mdr65KO-Wxac x w1118
## 2.120526e-14
## Mdr65KO x Excel_DF-Mdr65KO x w1118
## 4.396483e-14
## Mdr65KO x Bloom_DF-Mdr65KO x w1118
## 2.220446e-14
## Mdr65KO x Mdr65KO-Mdr65KO x w1118
## 2.087219e-14
## Mdr65KO x Bloom_DF-Mdr65KO x Excel_DF
## 5.803212e-02
## Mdr65KO x Mdr65KO-Mdr65KO x Excel_DF
## 1.063227e-02
## Mdr65KO x Mdr65KO-Mdr65KO x Bloom_DF
## 9.287183e-01Figure 5
And last but not least we have figure 5. This figure describes the quantification of the amount of insecticide (nitenpyram) in Mdr65KO flies compared to their control.
The data can be imported into R using the following code.
head(Fig5_HPLC_MS_Nitenpyram)## Rep Tissue Genotype Nitenpyram
## 1 1 Body Wxac 537186.33
## 2 2 Body Wxac 350397.70
## 3 3 Body Wxac 364973.92
## 4 1 Head Wxac 29739.99
## 5 2 Head Wxac 10007.03
## 6 3 Head Wxac 13127.84As before, the data needs for formatting before we can plot it.
options(scipen=999)
Tissue <- c("Body","Head")
annot <- data.frame(x = rep(-Inf, length(unique(Fig5_HPLC_MS_Nitenpyram$Tissue))), y = rep(Inf, length(unique(Fig5_HPLC_MS_Nitenpyram$Tissue))),
Tissue=as.character(unique(Fig5_HPLC_MS_Nitenpyram$Tissue)),labs=rev(LETTERS[1:length(unique(Fig5_HPLC_MS_Nitenpyram$Tissue))]))
Fig5_HPLC_MS_Nitenpyram <- Fig5_HPLC_MS_Nitenpyram[rev(order(Fig5_HPLC_MS_Nitenpyram$Genotype)),]
Fig5_HPLC_MS_Nitenpyram$Genotype <- factor(Fig5_HPLC_MS_Nitenpyram$Genotype,levels=c("Wxac","Mdr65KO"))
Fig5_HPLC_MS_Nitenpyram$Nitenpyram2 <- Fig5_HPLC_MS_Nitenpyram$Nitenpyram/10000Lastly we create the dotplot. There appears to be large differences between Mdr65KO and control. There are no errorbars because doing so with such a small sample size (n=3) is misleading.
gp <- ggplot(Fig5_HPLC_MS_Nitenpyram,aes(x=Genotype,y=Nitenpyram2,fill=Genotype))
gp <- gp+facet_grid(Tissue~.,scales="free_y")
gp <- gp+geom_dotplot(binaxis='y',stackdir="center",dotsize=2)
gp <- gp+scale_fill_manual(values=c("grey75","grey25"))
gp <- gp+ylab(label="Quantity of Nitenpyram\n(Area Under Chromatogram/10,000)\n")
gp <- gp+xlab("\nGenotype")
gp <- gp+geom_text(size=7,data=annot,aes(x, y, label=labs,fontface="bold",family="serif"), inherit.aes=FALSE,vjust=1.2,hjust=-.2)
gp <- gp+theme_bw(base_size=14,base_family="serif")
gp <- gp+theme(text=element_text(face="bold",family="serif"),panel.grid=element_blank(),axis.ticks.x=element_line(),
panel.border=element_rect(colour="black",fill=NA),strip.text=element_text(size=18,face="bold"),
strip.background=element_rect("grey95"),
axis.title=element_text(size=18),axis.text.x=element_text(angle=30,hjust=1),legend.position = "none")
print(gp)