###This script generates a plot of RATE-seq data for a gene of interest. ###It requires supplemental Tables S4 and S5 in the workspace #####Usage #Either a gene name or a systematic name is the sole argument to the function. #e.g. #1 #curve.generator("CTK1") # #e.g. #2 #curve.generator("YBR084C-A") #Created by B. Neymotin on May 20, 2014 #queries: TableS4 <- read.delim("TableS4.xls") TableS5 <- read.delim("TableS5.xls") mat.all <- merge(TableS4,TableS5) #put two matrices together curve.generator <- function(Syst){ time <- c(rep(3,3),rep(5,6),rep(7,6),rep(9,3),rep(11,3),rep(13,6),rep(20,3),rep(25,5),rep(30,2),rep(100,3)) #timepoints row.Gene <- which(mat.all[,42]==Syst) #the row with the gene name row.Syst <- which(mat.all[,1]==Syst) #the row with the systematic name row.Syst <- c(row.Gene,row.Syst) #combine the two. row.Syst <- row.Syst[1] if(is.na(row.Syst)){ x=paste(Syst,"is not included in the RATE-seq dataset") print(x) }else{ y.vals <- mat.all[row.Syst,2:41] #the yvals alpha.low <- mat.all[row.Syst,43] #the lower CI value for alpha alpha.real <- mat.all[row.Syst,44] #alpha alpha.high <- mat.all[row.Syst,45] #the higher CI value for alpha thalf.low <- signif(mat.all[row.Syst,49],digits=3) #the lower CI value for half-life thalf <- signif(mat.all[row.Syst,50],digits=3) #the value for half-life thalf.high <- signif(mat.all[row.Syst,51],digits=3) #the higher CI value for half-life Yeq <- mat.all[row.Syst,52] #the equilibrium value x.vals <- seq(0,100,1) y.vals.alpha <- Yeq*(1-exp(-alpha.real*(x.vals-2))) #make a curve with alpha y.vals.alpha.low <- Yeq*(1-exp(-alpha.low*(x.vals-2))) #make a curve with alpha.low y.vals.alpha.high <- Yeq*(1-exp(-alpha.high*(x.vals-2))) #make a curve with alpha.high #Plot data plot(time,y.vals,xlim=c(0,100),ylim=c(0,1.2*max(y.vals)),pch=17,col="purple",xlab="time (minutes)",ylab="counts (adjusted)",main=Syst) lines(x.vals,y.vals.alpha,col="purple") lines(x.vals,y.vals.alpha.low,col="black",lty=2) lines(x.vals,y.vals.alpha.high,col="black",lty=2) legend.name.alpha <- paste("alpha=",alpha.real,"min^-1","(95%CI ",alpha.low,"<->",alpha.high,"min^-1)") legend.name.thalf <- paste("half life=",thalf,"min","(95%CI ",thalf.low,"<->",thalf.high,"min)") legend("bottomright",c(as.expression(legend.name.alpha),as.expression(legend.name.thalf)),cex=0.8) } } ###Example plot: curve.generator("CTK1")