2012-05-01 80 views
2

这是我的数据和情节多个情节安排在基地情节中的R

nmar <- seq (1, 100, 5) 
position= rep(nmar, 5) 
n = length (nmar) 
chr = rep(1:5, each = n) 

mapdata <- data.frame (chr, position, 
snpname = paste("SNP-", 1:length (position), sep = "")) 
mapdata 


chr.lab = 1 ; mbar.col = "blue" 
layout(matrix(c(1,1,2),nc=1)) # works for two but I need to extend it to 
     n (which is level of chr = 5) 

# plot level 1 
mapdata1 <- mapdata[mapdata$chr == 1,] 
m <- par()$mar 
m[1] <- m[3] <- 0 
par(mar=m) 
# Set the limits of the plot 
plot(mapdata1$position,mapdata1$position-mapdata1$position, type="n", 
    axes=FALSE, 
xlab="", ylab="Chromsome", yaxt="n") 

polygon(
    c(0,max(mapdata1$position + 0.08 * max(mapdata1$position)),max(mapdata1$position)+ 
    0.08 * max(mapdata1$position),0), 
    .2*c(-0.3,-0.3,0.3,0.3), 
    col= mbar.col 
) 
segments(mapdata1$position, -.3, mapdata1$position, .3) 
text(mapdata1$position, -.7, mapdata1$snpname, srt = 90, cex.lab = chr.lab) 
text(mapdata1$position, .7, mapdata1$position,cex.lab = chr.lab) 
text(0, labels = c("Chr 2")) 

第二级

# plot level 2 
mapdata2 <- mapdata[mapdata$chr == 2,] 
m <- par()$mar 
m[1] <- m[3] <- 0 
par(mar=m) 
# Set the limits of the plot 
plot(mapdata2$position,mapdata2$position-mapdata1$position, type="n", axes=FALSE, 
xlab="", ylab="Chromsome", yaxt="n") 

polygon(
    c(0,max(mapdata2$position + 0.08 * max(mapdata2$position)),max(mapdata2$position)+ 
0.08 * max(mapdata2$position),0), 
    .2*c(-0.3,-0.3,0.3,0.3), 
    col= mbar.col 
) 
segments(mapdata2$position, -.3, mapdata2$position, .3) 
text(mapdata2$position, -.7, mapdata2$snpname, srt = 90, cex.lab = chr.lab) 
text(mapdata2$position, .7, mapdata2$position,cex.lab = chr.lab) 
text(0, labels = c("Chr 2")) 

输出 enter image description here

(1)如何能我自动化了n个级别的过程 - 将类似的绘图扩展到n个级别的chr (2)如果您看到具有相同规格的barsize具有c绞死,可能是由于不同的地块面积。我怎样才能调整它,使所有的地块相同?

+0

有什么特别的理由想坚持到基地地块,这样的函数取值CHR的相反,你可以修改吗? ggplot2和lattice为创建这些类型的复合图解提供了简单的解决方案,请参阅这些链接以获得大量示例http://learnr.wordpress.com/2009/06/28/ggplot2-version-of-figures-in-lattice-multivariate -data-visualization-with-r-part-1/http://learnr.wordpress.com/2009/06/29/ggplot2-version-of-figures-in-lattice-multivariate-data-visualization-with-r -part-2/http://learnr.wordpress.com/2009/07/02/ggplot2-version-of-figures-in-lattice-multivariate-data-visualization-with-r-part-3/ –

+0

我同意ggplot2和格子更好,但我觉得他们很难操纵,不是吗? – SHRram

+0

我发现ggplot2很容易使用,一旦你掌握了它。 Imo它确实值得它切换到ggplot2,或者说格子,虽然在使用它们两个之后我更喜欢ggplot2。 –

回答

2

ggplot绝对是这里的路。但如果你真的想坚持基地plot,那么这个功能会工作:

plot.as.stack= function(mapdata1, mbar.col = "blue"){ 
    # mapdata1 <- mapdata[mapdata$chr == chr,] 
    m <- par()$mar 
    m[1] <- m[3] <- 0 
    par(mar=m) 
    # Set the limits of the plot 
    plot(mapdata1$position,mapdata1$position-mapdata1$position, type="n", 
     axes=FALSE, 
    xlab="", ylab="Chromsome", yaxt="n") 

    polygon(
     c(0,max(mapdata1$position + 0.08 * max(mapdata1$position)),max(mapdata1$position)+ 
     0.08 * max(mapdata1$position),0), 
     .2*c(-0.3,-0.3,0.3,0.3), 
     col= mbar.col 
    ) 
    segments(mapdata1$position, -.3, mapdata1$position, .3) 
    text(mapdata1$position, -.7, mapdata1$snpname, srt = 90, cex.lab = chr.lab) 
    text(mapdata1$position, .7, mapdata1$position,cex.lab = chr.lab) 
    text(0, labels = paste("Chr",unique(mapdata1$chr))) 
} 

# Example Run. 
par(mfrow=c(length(unique(mapdata$chr)),1)) 
x=by(mapdata,factor(mapdata$chr),plot.as.stack) # Assigned to x to prevent output 
par(mfrow=c(1,1)) 

正如你所看到的,我只是把你的代码,把它变成一个函数,然后在它跑by。请注意,这将运行的所有上的功能chr的级别。

plot.as.stack = function(chr){ 
    mapdata1 <- mapdata[mapdata$chr == chr,] 
    ... 
} 

然后用CHR值运行函数:

par(mfrow=c(5,1)) 
sapply(1:5,plot.as.stack) 
par(mfrow=c(1,1)) 
+0

我迫不及待地接受你的答案,谢谢! – SHRram