Mercurial > repos > xuebing > sharplabtool
comparison tools/mytools/cdf-old-not-used/cdf.py @ 0:9071e359b9a3
Uploaded
author | xuebing |
---|---|
date | Fri, 09 Mar 2012 19:37:19 -0500 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:9071e359b9a3 |
---|---|
1 ''' | |
2 python wrapper for cdf plot | |
3 ''' | |
4 | |
5 import os,sys | |
6 | |
7 filename = sys.argv[1] | |
8 output = sys.argv[2] | |
9 c1, c2 = sys.argv[3].strip().split(',') | |
10 label1,label2 = sys.argv[5].strip().split(',') | |
11 title = sys.argv[6].strip() | |
12 log = sys.argv[4] | |
13 | |
14 rf = open('tmp.r','w') | |
15 | |
16 rf.write("options(warn=-1)\n") | |
17 | |
18 rf.write("x <- read.table('"+filename+"')\n") | |
19 | |
20 rf.write("mycdf <- function(x1,x2) {\n") | |
21 rf.write(" sx1 <- sort(x1)\n") | |
22 rf.write(" sy1 <- c(1:length(sx1))/length(sx1)\n") | |
23 rf.write(" sx2 <- sort(x2)\n") | |
24 rf.write(" sy2 <- c(1:length(sx2))/length(sx2)\n") | |
25 rf.write(" mi <- min(c(x1,x2))\n") | |
26 rf.write(" ma <- max(c(x1,x2))\n") | |
27 rf.write(" sx1 <- c(mi,sx1,ma)\n") | |
28 rf.write(" sx2 <- c(mi,sx2,ma)\n") | |
29 rf.write(" sy1 <- c(0,sy1,1)\n") | |
30 rf.write(" sy2 <- c(0,sy2,1)\n") | |
31 rf.write(" pv <- ks.test(x1,x2)$p.value\n") | |
32 rf.write(" list(sx1,sy1,sx2,sy2,pv)\n") | |
33 rf.write("}\n") | |
34 | |
35 rf.write("list1 <- x[,"+c1+"]\n") | |
36 rf.write("list2 <- x[,"+c2+"]\n") | |
37 | |
38 if log == 'log': | |
39 rf.write("list1 <- log2(1+list1) \n") | |
40 rf.write("list2 <- log2(1+list2) \n") | |
41 #rf.write("list1[list1 == -Inf] <- min(list1[list1 != -Inf])\n") | |
42 rf.write("res <- mycdf(list1,list2)\n") | |
43 | |
44 rf.write("pdf('"+output+"')\n") | |
45 rf.write("plot(res[[1]],res[[2]],type='l',lty=1,lwd=2,col='red',main=paste('"+title+"',format(res[[5]], digit=2,scientific = T),sep=' , p='),xlab='value',ylab='cumulative frequency')\n") | |
46 rf.write("lines(res[[3]],res[[4]],type='l',lty=2,lwd=2,col='blue')\n") | |
47 rf.write("label1 <- paste('"+label1+"','( n1=',length(list1),')')\n") | |
48 rf.write("label2 <- paste('"+label2+"','( n2=',length(list2),')')\n") | |
49 rf.write("legend('bottomright',c(label1,label2),col=c('red','blue'), lty=1:2, lwd=2, bty='n')\n") | |
50 rf.write("boxplot(list1,list2,names=c(label1,label2))\n") | |
51 rf.write("dev.off()\n") | |
52 | |
53 rf.close() | |
54 os.system("R --vanilla < tmp.r") | |
55 os.system('rm tmp.r') |