Question: adding marker(snp-position) to the regional association plot
0
gravatar for muralinmars
4.4 years ago by
muralinmars70
France
muralinmars70 wrote:

Hi all,

I have a dataset for siginificant cpg probes (with p-value) associated with snps in a gene region... I managed to plot the data as a regional plot  (in a span of 1mb region)..

Now I would like to add marker or SNP position to the plot... something like this.. (figure A)

http://www.cell.com/cms/attachment/2016737868/2037292392/gr2.jpg

here is the r code(from another source) ... that takes in the pvalues and the gene positions..

rectarrows <- function(x0,y0,x1,y1,height,length,...) {
  lwd=par("lwd")
  l0=height*(y1-y0)/sqrt((x1-x0)^2+(y1-y0)^2)
  l1=height*(x1-x0)/sqrt((x1-x0)^2+(y1-y0)^2)
  d0=length*(y1-y0)/sqrt((x1-x0)^2+(y1-y0)^2)
  d1=length*(x1-x0)/sqrt((x1-x0)^2+(y1-y0)^2)
  polygon(x=c(x0+l0,x1+l0-d1,x1,x1-l0-d1,x0-l0),y=c(y0-l1,y1-l1-d0,y1,y1+l1-d0,y0+l1),...)
}

pval<-runif(6000, 0, 1)
logPval<--log(pval,base=10)
pos=1:6000
rsID<-paste("rs",1:6000,sep="")
data<-as.data.frame(cbind(pos,rsID,pval,logPval))

genes <- c( "NFIA", "KANK4", "DOCK7", "L1TD1")
gene.start<- c(4, 350, 1000, 2100)
gene.end<- c(3004, 3350, 4000, 5010)
genes.pos<- as.data.frame(cbind(genes, gene.start, gene.end))

genes.pos[,1]=as.character(genes.pos[,1])
genes.pos[,2]=as.numeric(as.character(genes.pos[,2]))
genes.pos[,3]=as.numeric(as.character(genes.pos[,3]))

lay=layout(matrix(seq(2),2,1,byrow=TRUE),heights=c(2000,1000))
par(mar=c(0,4.1,4.1,2.1))
plot(x=as.numeric(data$pos),y=as.numeric(as.character(data$logPval)),xaxt="n")
par(mar=c(5.1,4.1,0.5,2.1))
plot(0,0,type="n",xlim=c(min(pos),max(pos)),ylim=c(0,length(genes)*1.1),yaxt="n")
lapply(seq(length(genes)), function(i) {
  rectarrows(genes.pos$gene.start[i],i,genes.pos$gene.end[i],i,col="blue",height=0.1,length=100) # change to length=1000 or length=10000 for large x-

axes, as this is the length of the 'arrow' part of the rectangle
  text(mean(c(genes.pos$gene.start[i],genes.pos$gene.end[i])),i+0.3,genes.pos$genes[i],cex=0.7)
})

 

NOW HOW CAN EDIT THIS TO HIGHLIGHT THE SNP MARKER IN THE PLOT (AS SHOWN IN FIGUREa)

thanks

M

plot snp sequencing R • 1.3k views
ADD COMMENTlink written 4.4 years ago by muralinmars70
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 713 users visited in the last hour