This is what I'm looking for - Ill have a go! Thanks
Hi Folks
I wanted to make a plot in R and Im struggling with it!
I was hoping to use a BED file to plot blocks of equal height but varying length (peak positions) along the x axis (genomic position) the reason I want equal height is so I could plot using color intensity based on the fold enrichment values from the MACs peak caller. Exactly like the UCSC genome browser "dense" view. Can somebody give me an example of this??
Many thanks in advance
Basically like this track on hg19
2 answers
This is a pretty low level way of achieving it, if I correctly understand:
## Dummy data in ~bed format
n<- 10
chrom<- rep('chr1', n)
start<- seq(1, 1000, length.out= n)
end<- start + rpois(n, 50)
logfc<- rgamma(n, 2)
## Prepare colours based on logFC
pal<- colorRampPalette(c('blue', 'red'))
cols<- paste0(pal(n)[as.numeric(cut(logfc, breaks = n))], 99)
## Plot
plot(x= start, y= rep(0, n), type= 'n', bty= 'n', yaxt= 'n',
ylab= '', xlab= 'Position', xlim= range(c(start, end)), ylim= c(0, 1))
rect(xleft= start, xright= end, ybottom= 0, ytop= 0.1, border= NA, col= cols)
text(x= start, y= 0.12, labels= sprintf("%.2f", logfc), adj= c(0,0))
## Legend
nlg<- 5
lgcols<- paste0(pal(nlg)[1:nlg], 99)
lg<- round(seq(min(logfc), max(logfc), length.out= nlg), 1)
legend('topleft', pch= 15, bty= 'n', col= lgcols, legend= lg, pt.cex= 2, cex= 1)

Ok so I got it to work with my data and it looks good - thank you - However I was hoping to get a color legend embedded also. Do you have any Idea how I would do that? I have looked at plot3D
See edit with legend. it's not very nice because the numbers are not nice & round, but maybe it gives you a hint...
Maybe just export the UCSC browser shot directly from an existing session, or use a fetch tool, if automation is needed.
We are using a custom genome (not on UCSC) and as far as i know we cant add a custom genome to it. We like our previous outputs using R which is why we were hoping to use it again
Log in to answer this question.