-
Notifications
You must be signed in to change notification settings - Fork 1
/
metaplot.R
83 lines (66 loc) · 3.01 KB
/
metaplot.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
#metaplot
args = commandArgs(trailingOnly=TRUE)
samplelist=read.table("samples_list.txt",sep='\t')
#n_rain=nrow(samplelist)
#cl <- rainbow(n_rain)
metafile1=read.table(paste(samplelist[1,1],"_CG.matrix.gz",sep=''),skip=1)
metafile2=read.table(paste(samplelist[1,1],"_CHG.matrix.gz",sep=''),skip=1)
metafile3=read.table(paste(samplelist[1,1],"_CHH.matrix.gz",sep=''),skip=1)
y_cg=max(apply(as.matrix(metafile1[,7:86]),2,mean,na.rm=T))*1.2*100
y_chg=max(apply(as.matrix(metafile2[,7:86]),2,mean,na.rm=T))*1.2*100
y_chh=max(apply(as.matrix(metafile3[,7:86]),2,mean,na.rm=T))*1.2*100
outfile=paste("./metaplot_CG.png",sep="")
bitmap(outfile,type="png16m",height=3.5,width=3.5,res=450)
plot(NA,ylab='CG methylation levels (%)',xlab=NA,type='l',cex.main = 1,xlim=c(0,ncol(metafile1)),ylim=c(0,y_cg),xaxt="n",cex=1.5)
#palette()[1] "black" "red" "green3" "blue" "cyan" "magenta" "yellow" "gray"
for (i in 1:nrow(samplelist)){
name=samplelist[i,1]
file=paste(name,"_CG.matrix.gz",sep='')
metafile=read.table(file,sep='\t',skip=1)
av1=apply(as.matrix(metafile[,7:86]),2,mean,na.rm=T)
lines(av1*100,xlab=NA,type='l',cex.main = 1,xaxt="n",cex=1.5,col=palette()[i])
}
xx<-c(20,40,60)
ind<-c("TSS",args[1],"TES")
abline(v=20,lty=2,col='grey70')
abline(v=60,lty=2,col='grey70')
axis(1,at=xx,labels=ind)
legend("topright", legend = t(samplelist[1]), col = t(palette()), lwd = 1,cex = 0.5)
dev.off()
outfile=paste("./metaplot_CHG.png",sep="")
bitmap(outfile,type="png16m",height=3.5,width=3.5,res=450)
plot(NA,ylab='CHG methylation levels (%)',xlab=NA,type='l',cex.main = 1,xlim=c(0,ncol(metafile2)),ylim=c(0,y_chg),xaxt="n",cex=1.5)
#palette()[1] "black" "red" "green3" "blue" "cyan" "magenta" "yellow" "gray"
for (i in 1:nrow(samplelist)){
name=samplelist[i,1]
file=paste(name,"_CHG.matrix.gz",sep='')
metafile=read.table(file,sep='\t',skip=1)
av1=apply(as.matrix(metafile[,7:86]),2,mean,na.rm=T)
lines(av1*100,xlab=NA,type='l',cex.main = 1,xaxt="n",cex=1.5,col=palette()[i])
}
xx<-c(20,40,60)
ind<-c("TSS",args[1],"TES")
abline(v=20,lty=2,col='grey70')
abline(v=60,lty=2,col='grey70')
axis(1,at=xx,labels=ind)
legend("topright", legend = t(samplelist[1]), col = t(palette()), lwd = 1,cex = 0.5)
dev.off()
#CHH
outfile=paste("./metaplot_CHH.png",sep="")
bitmap(outfile,type="png16m",height=3.5,width=3.5,res=450)
plot(NA,ylab='CHH methylation levels (%)',xlab=NA,type='l',cex.main = 1,xlim=c(0,ncol(metafile3)),ylim=c(0,y_chh),xaxt="n",cex=1.5)
#palette()[1] "black" "red" "green3" "blue" "cyan" "magenta" "yellow" "gray"
for (i in 1:nrow(samplelist)){
name=samplelist[i,1]
file=paste(name,"_CHH.matrix.gz",sep='')
metafile=read.table(file,sep='\t',skip=1)
av1=apply(as.matrix(metafile[,7:86]),2,mean,na.rm=T)
lines(av1*100,xlab=NA,type='l',cex.main = 1,xaxt="n",cex=1.5,col=palette()[i])
}
xx<-c(20,40,60)
ind<-c("TSS",args[1],"TES")
abline(v=20,lty=2,col='grey70')
abline(v=60,lty=2,col='grey70')
axis(1,at=xx,labels=ind)
legend("topright", legend = t(samplelist[1]), col = t(palette()), lwd = 1,cex = 0.5)
dev.off()