求助,用R语言plot函数做manhattan散点图

COS论坛 | 统计之都 COS论坛 | 统计之都 软件应用 S-Plus & R语言 求助,用R语言plot函数做manhattan散点图

该主题包含 2 条回复,3个帖子,最后由  dustykl3 周 之前 更新。

查看 3 个帖子 - 1 到 3(总计 3 个)
  • 作者
    帖子
  • 1 楼

    qingcanfan
    Participant

    我根据网上的程序画出manhattan散点图,命令为

    manhattan = function(dataframe, colors=c(“black”,”#666666″,”#CC6600″), pch=20,ymax=”max”, cex.x.axis=1, limitchromosomes=1:23,suggestiveline=-log10(3.0E-08),genomewideline=-log10(1.8-06), annotate=NULL, …) {

    d=dataframe
    if (!(“CHR” %in% names(d) & “BP” %in% names(d) & “P” %in% names(d))) stop(“Make sure your data frame contains columns CHR, BP, and P”)

    if (any(limitchromosomes)) d=d[d$CHR %in% limitchromosomes, ]
    d=subset(na.omit(d[order(d$CHR, d$BP), ]), (P>0 & P<=1)) # remove na’s, sort, and keep only 0<P<=1
    d$logp = -log10(d$P)
    d$pos=NA
    ticks=NULL
    lastbase=0
    colors <- rep(colors,max(d$CHR))[1:max(d$CHR)]
    if (ymax==”max”) ymax<-ceiling(max(d$logp))
    if (ymax<8) ymax<-8

    numchroms=length(unique(d$CHR))
    if (numchroms==1) {
    d$pos=d$BP+1
    ticks=floor(length(d$pos))/2+1
    } else {
    for (i in unique(d$CHR)) {
    if (i==1) {
    d[d$CHR==i, ]$pos=d[d$CHR==i, ]$BP
    } else {
    lastbase=lastbase+tail(subset(d,CHR==i-1)$BP, 1)
    d[d$CHR==i, ]$pos=d[d$CHR==i, ]$BP+lastbase
    }
    ticks=c(ticks, d[d$CHR==i, ]$pos[floor(length(d[d$CHR==i, ]$pos)/2)+1])
    }
    }

    if (numchroms==1) {
    with(d, plot(pos, logp, ylim=c(0,ymax), ylab=expression(-log[10](italic(p))), xlab=paste(“Chromosome”,unique(d$CHR),”position”), pch=20))
    } else {
    with(d, plot(pos, logp, ylim=c(0,ymax), ylab=expression(-log[10](italic(p))), xlab=”Chromosome”, xaxt=”n”, type=”n”,pch=20))
    axis(1, at=ticks, lab=unique(d$CHR), cex.axis=cex.x.axis)
    icol=1
    for (i in unique(d$CHR)) {
    with(d[d$CHR==i, ],points(pos, logp, col=colors[icol], pch=20))
    icol=icol+1
    }
    }

    if (!is.null(annotate)) {
    d.annotate=d[which(d$SNP %in% annotate), ]
    with(d.annotate, points(pos, logp, col=”green3″, pch=20))
    }

    if (genomewideline) abline(h=genomewideline, col=”red”)
    if (potentialline) abline(h=potentialline, col=”green”)
    }

    现在图已经汇出,但是我想在x轴每个数值间加一个空白间隔,如下图跪求各位大神帮忙解决。
    http://www.google.com.hk/imgres?imgurl=http://www.nature.com/ng/journal/v41/n9/images/ng.429-F1.jpg&imgrefurl=http://www.nature.com/ng/journal/v41/n9/fig_tab/ng.429_F1.html&h=403&w=900&sz=87&tbnid=fBPZNz9Y2JRD6M:&tbnh=62&tbnw=139&prev=/search%3Fq%3Dgwas%2Bmanhattan%25E5%259B%25BE%26tbm%3Disch%26tbo%3Du&zoom=1&q=gwas+manhattan%E5%9B%BE&usg=__DCE8nvzwcCGZMHLLswmKre6oHuA=&docid=-XvbA2hEiRrCJM&hl=zh-CN&sa=X&ei=oXM9UfOdKo2jiget_YCACQ&ved=0CEEQ9QEwBQ&dur=876

    2 楼

    jianggl2000
    Participant

    lastbase=lastbase+tail(subset(d,CHR==i-1)$BP, 1) * 1.5

    如果间隙太大,修正1.5的大小
    或者去除*1.5,换为“+一个固定值”

    3 楼

    dustykl
    Participant

    回复 2 楼jianggl2000
    你好,请问加空白间隙这个解决了吗?我在做图的过程中也没有找到六空白间隙的方法。

查看 3 个帖子 - 1 到 3(总计 3 个)

您必须先登录才能回复该主题。