========================================================
MF batch: XBATCHX

library(knitr)
library(dplyr)
library(ggplot2)
library(gridExtra)
library(ggthemes)
data<-rio::import('data.csv')
models<-lapply(split(data,data$dye),function(x){nls(counts ~ SSfpl(pH, Bottom, Top, pKa, Slope), x)})
model.names<-models  %>% {names(.) %>% split(.,.)}
pHvals <- seq(min(data$pH), max(data$pH), length.out=150) 
linedata<-model.names %>% 
  lapply(.,{.  %>%  
      models[[.]]  %>% 
      predict(.,data.frame(pH=pHvals)) %>% 
      .[1:150]})
titrlines<-lapply(model.names,{.  %>% data.frame(val=linedata[[.]],dye=.,pH=pHvals)}) %>% bind_rows()
cols <- c("CL" = "cornflowerblue","PR" = "red")

titrPlot<-ggplot(data,aes(x=pH,y=counts,colour=dye,fill=dye))+geom_jitter(aes(fill=dye),size=2,shape=16,alpha=.7)+
geom_line(data=titrlines,aes(x=pH,y=val,colour=dye))+
theme_few()+ggtitle('Titration Curve')+
guides(colour = FALSE,fill=FALSE) +
scale_colour_manual(values = cols)


## lines for effective range
min10<-function(z){min(z) + ((max(z) - min(z))/10)}
max10<-function(z){max(z) - ((max(z) - min(z))/10)}
LowLine<-function(u){b<-min10(u);which(abs(u - b) == min(abs(u - b)))}
HighLine<-function(u){b<-max10(u);which(abs(u - b) == min(abs(u - b)))}
#
pH.index<-c(lapply(linedata,LowLine),lapply(linedata,HighLine)) %>% Reduce('rbind',x=.)

EfRange<-data.frame(vline=pHvals[pH.index],dye=factor(c('CL','PR')))

 gainData<-models %>% 
   names(.) %>% 
   split(.,.) %>% 
    lapply(.,{. %>%
  data.frame(gain=diff(linedata[[.]][1:150])/diff(pHvals)/1000,pH=pHvals[-1],dye=.)}) %>% 
    bind_rows() %>%
    mutate(.data=.,dye=factor(dye))


gain_plot<-ggplot(gainData,aes(x=pH,y=gain,colour=dye)) + geom_line(size=2) +theme_few() +ggtitle('Gain Curve') + geom_vline(data=EfRange,aes(xintercept=vline,colour=dye),linetype = 'dotted')+scale_colour_manual(values = cols)
grid.arrange(titrPlot,gain_plot, nrow=1, ncol=2,top='pKA')

Effective range

EfRange$stat<-c(rep('Low.pH',2),rep('High.pH',2))
fix<-c('CL'='Clear','PR'='Phenol Red')
EfRange$dye<-fix[EfRange$dye]

r kable(t(xtabs(vline~.,data=EfRange)),format='markdown')

Clear pKa information

Clear.V50<-cbind(signif(summary(models[['CL']])$parameters,digits=5),
signif(confint(models[['CL']]),digits=4))

r kable(Clear.V50,format='markdown')

Phenol Red pKa information

Phenol.Red.V50<-cbind(signif(summary(models[['PR']])$parameters,digits=5),
signif(confint(models[['PR']]),digits=4))

r kable(Phenol.Red.V50,format='markdown')



JARS3N/PipeFish documentation built on May 7, 2019, 6:47 a.m.