library(preprocessCore)


results.thread1 <- NULL
results.thread2 <- NULL
results.thread3 <- NULL
for (ncol in c(10,20,50,100,200)){
  for (nrow in c(10000,100000,1000000,1000000)){
      X <- matrix(rexp(nrow*ncol),nrow,ncol)
      Sys.setenv(R_THREADS=1)	
      results.thread1 <- c(results.thread1,list(system.time(Y <- normalize.quantiles(X))))	
      Sys.setenv(R_THREADS=2)
      results.thread2 <- c(results.thread2,list(system.time(Y <- normalize.quantiles(X))))
      Sys.setenv(R_THREADS=3)
      results.thread3 <- c(results.thread3,list(system.time(Y <- normalize.quantiles(X))))
   }
}


results.thread1 <- do.call("rbind",results.thread1)	
results.thread2 <- do.call("rbind",results.thread2)	
results.thread3 <- do.call("rbind",results.thread3)	

ncol <- rep(c(10,20,50,100,200),c(4,4,4,4,4))
nrow <- rep(c(10000,100000,1000000,1000000),5)


png("preprocessCore_QN_times_%03d.png",height=1024,width=1024)
plot(ncol[nrow==10000],results.thread1[nrow==10000,3],log="xy",type="l",ylab="Time in seconds",xlab="Number of Columns",main="Quantile Normalization Run times (10,000 rows)")
lines(ncol[nrow==10000],results.thread2[nrow==10000,3],col="blue")
lines(ncol[nrow==10000],results.thread3[nrow==10000,3],col="red")
legend(10,1.0,c("R_THREADS=1","R_THREADS=2","R_THREADS=3"),
	col=c("black","blue","red"),lty=1,lwd=2)


plot(ncol[nrow==100000],results.thread1[nrow==100000,3],log="xy",type="l",ylab="Time in seconds",xlab="Number of Columns",main="Quantile Normalization Run times (100,000 rows)",ylim=c(0.3,12))
lines(ncol[nrow==100000],results.thread2[nrow==100000,3],col="blue")
lines(ncol[nrow==100000],results.thread3[nrow==100000,3],col="red")
legend(10,10,c("R_THREADS=1","R_THREADS=2","R_THREADS=3"),
	col=c("black","blue","red"),lty=1,lwd=2)


plot(ncol[nrow==1000000],results.thread1[nrow==1000000,3],log="xy",type="l",ylab="Time in seconds",xlab="Number of Columns",main="Quantile Normalization Run times (1,000,000 rows)",ylim=c(4,150))
lines(ncol[nrow==1000000][c(1,3,5,7,9)],(results.thread2[nrow==1000000,3][c(1,3,5,7,9)] + results.thread2[nrow==1000000,3][c(2,4,6,8,10)])/2,col="blue")
lines(ncol[nrow==1000000][c(1,3,5,7,9)],(results.thread3[nrow==1000000,3][c(1,3,5,7,9)] + results.thread3[nrow==1000000,3][c(2,4,6,8,10)])/2,col="red")
legend(10,100,c("R_THREADS=1","R_THREADS=2","R_THREADS=3"),
	col=c("black","blue","red"),lty=1,lwd=2)


plot(ncol[nrow==100000],results.thread2[nrow==10000,3]/results.thread1[nrow==10000,3],type="l",ylim=c(0.4,1),lty=1,col="brown",ylab="Ratio of 2 Thread to Single Thread run time",xlab="Number of Columns",main="Comparing RTHREADS=2 to RTHREADS=1")
lines(ncol[nrow==100000],results.thread2[nrow==100000,3]/results.thread1[nrow==100000,3],lty=1,col="orange")
lines(ncol[nrow==100000],(results.thread2[nrow==1000000,3]/results.thread1[nrow==1000000,3])[c(2,4,6,8,10)],lty=1,col="green")
legend(125,1,c("rows=10,000","rows=100,000","rows=1,000,000"),lty=1,col=c("brown","orange","green"))


plot(ncol[nrow==100000],results.thread3[nrow==10000,3]/results.thread1[nrow==10000,3],type="l",ylim=c(0.4,1),lty=1,col="brown",ylab="Ratio of 3 Thread to Single Thread run time",xlab="Number of Columns",main="Comparing RTHREADS=3 to RTHREADS=1")
lines(ncol[nrow==100000],results.thread3[nrow==100000,3]/results.thread1[nrow==100000,3],lty=1,col="orange")
lines(ncol[nrow==100000],(results.thread3[nrow==1000000,3]/results.thread1[nrow==1000000,3])[c(2,4,6,8,10)],lty=1,col="green")
legend(125,1,c("rows=10,000","rows=100,000","rows=1,000,000"),lty=1,col=c("brown","orange","green"))


plot(ncol[nrow==100000],results.thread3[nrow==10000,3]/results.thread2[nrow==10000,3],type="l",ylim=c(0.4,1.5),lty=1,col="brown",ylab="Ratio of 3 Thread to 2 Thread run time",xlab="Number of Columns",main="Comparing RTHREADS=3 to RTHREADS=2")
lines(ncol[nrow==100000],results.thread3[nrow==100000,3]/results.thread2[nrow==100000,3],lty=1,col="orange")
lines(ncol[nrow==100000],(results.thread3[nrow==1000000,3]/results.thread2[nrow==1000000,3])[c(2,4,6,8,10)],lty=1,col="green")
legend(125,1.5,c("rows=10,000","rows=100,000","rows=1,000,000"),lty=1,col=c("brown","orange","green"))
dev.off()


sink("Session")
sessionInfo()
sink()


save(results.thread1, results.thread2, results.thread3,file="Runtimes.Rda")
q(save="no")