MotifStack: motif 可视化
最近大量跑chip-seq,看到一篇2016Cell的文章《Cistrome and Epicistrome Features Shape the Regulatory DNA Landscape》感觉图3好惊艳。
找了下材料和方法,发现是一个叫MotifStack(http://www.bioconductor.org/packages/release/bioc/html/motifStack.html)的包画的。
什么是motif?
Motif是在生物学中是一个基于数据的数学统计模型,典型的是一段sequence也可以是一个结构。如转录因子倾向于结合某些特定的序列。把这个R包和数据下载下来自己画一下,它motif的测试数据是pcm格式的文本。MotifStack 还有一些motif格式转换的函数。
motif的格式是什么?
Chip-seq鉴定出的motif,是一个ATCG的序列矩阵。对于这个字母矩阵,目前主流Motif的序列格式主要有JASPAR\MEME\RAW PFM。此外JASPAR数据库 (http://jaspar.genereg.net/) 提供了转录因子与DNA结合位点motif最全面的公开数据,共收集了脊椎动物、植物、昆虫、线虫、真菌和尾索动物六大类不同类生物的数据。
1 Single motif log
pcm <- read.table("/Users/zt/learn/软件包/motifStack/test_data/bin_SOLEXA.pcm",header = F)# 创造Motif对象pcm <- pcm[,3:ncol(pcm)]rownames(pcm) <- c("A","C","G","T")motif <- new("pcm", mat=as.matrix(pcm), name="bin_SOLEXA")opar<-par(mfrow=c(4,1))plot(motif)#plot the logo with same heightplot(motif, ic.scale=FALSE, ylab="probability")#try a different fontplot(motif, font="mono,Courier")#try a different font and a different color groupmotif@color <- colorset(colorScheme='basepairing')plot(motif,font="Times")
2 Affinity motif log
需要画出双链。
motif<-matrix( c( .846, .631, .593, .000, .000, .000, .434, .410, 1.00, .655, .284, .000, .000, .771, .640, .961, .625, .679, .773, 1.00, 1.00, .000, .573, .238, .397, 1.00, 1.00, .000, .298, 1.00, 1.00, .996, 1.00, 1.00, 1.00, .228, .000, 1.00, 1.00, .597, .622, .630, .000, 1.00, 1.00, .871, .617, 1.00, .701, .513, .658, .000, .000, .247, .542, 1.00, .718, .686, .000, .000, .000, .595, .437, .970 ), nrow=4, byrow = TRUE)rownames(motif) <- c("A", "C", "G", "T")motif<-new("psam", mat=motif, name="affinity logo")plot(motif)
3 Stacked motif
当要画多个motif以堆积起来,对不同的motif进行距离计算构建进化树。导入的数据文件是单个motif一个pcm文件,通过file.path命令读取整个目录
motifs<-importMatrix(dir(file.path("/Users/zt/learn/软件包/motifStack/test_data"),"pcm$", full.names = TRUE))motifStack(motifs, layout="stack", ncex=1.0)# 建立进化树motifStack(motifs, layout="tree")
4 Circle motif
设置layout的格式,可以设置成环形,有点像circles
motifStack(motifs, layout="radialPhylog", circle=0.5, cleaves = 0.3, clabel.leaves = 0.3, col.bg=rep(color, each=5), col.bg.alpha=0.3, col.leaves=rep(color, each=5), col.inner.label.circle=rep(color, each=5), inner.label.circle.width=0.05, col.outer.label.circle=rep(color, each=5), outer.label.circle.width=0.1, circle.motif=1.2, angle=350)
相关推荐:
- NCBI在线BLAST使用方法与结果详解 2938
- 神经网络术语:Epoch、Batch Size和迭代 527
- Consed的安装与使用教程 465
- 陈连福的NGS生物信息学培训教材V2.1 277
- WGCNA分析使用教程 272
最新创建圈子
-
原料药研发及国内外注册申报
2019-01-25 10:41圈主:caolianhui 帖子:33 -
制药工程交流
2019-01-25 10:40圈主:polysciences 帖子:30 -
健康管理
2019-01-25 10:40圈主:neuromics 帖子:20 -
发酵技术
2019-01-25 10:39圈主:fitzgerald 帖子:17 -
医学肿瘤学临床试验
2019-01-25 10:39圈主:bma 帖子:58