diff genomeview_notused @ 20:16ba480adf96

Uploaded

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/genomeview_notused	Sat Mar 31 08:31:22 2012 -0400
@@ -0,0 +1,45 @@
+
+caloffset = function(genome){
+    total_len = sum(as.numeric(genome[,2]))
+    offset = 0
+    for (i in 1:nrow(genome)) {
+        offset = c(offset,offset[i]+genome[i,2])        
+    }
+    offset    
+}
+
+coverage = function(intervals,genome,offset,resolution) {
+
+    nChr = length(offset) - 1
+    total_len = offset[nChr+1]
+    nbin = as.integer(total_len / resolution)
+
+	pos = numeric(0)
+    cov = numeric(0)#coverage
+    col = numeric(0)#color
+    for (i in 1:nChr) {
+        d = x[x[,1]==as.character(genome[i,1]),2:3]
+        d = ceiling(((d[,1]+d[,2])/2+offset[i])*nbin/total_len)
+        t = table(d)
+		pos = c(pos,as.numeric(row.names(t)))
+        cov = c(cov, as.numeric(t))
+        col = c(col,numeric(length(t))+i)
+    }
+    list(nbin=nbin,pos=pos,cov=cov,col=col)
+}
+
+# plot coverage
+# res = genomeView(x,genome,100000)
+plotcov = function(res,genome,offset,title,uselog) {
+	if (uselog == 'log'){
+		res$cov = log10(res$cov+1)
+	}
+    ymax = max(res$cov)
+	par(mar=c(5,5,5,1))
+    plot(res$pos,res$cov,type='h',cex=0.1,cex.axis=2,cex.lab=2,cex.main=3,col=res$col,xaxt='n',main=title,xlab='chromosome',ylab='coverage',frame.plot=F,xlim=c(0,res$nbin),ylim=c(0,ymax))
+    xticks = numeric(nrow(genome))
+    for (i in 1:nrow(genome)){
+       xticks[i] = (offset[i]+offset[i+1])/2*res$nbin/offset[length(offset)]
+    }
+    mtext(genome[,1],side=1,at=xticks,adj=1,las=2,col=seq(nrow(genome)))
+}