| 查看: 1099 | 回复: 0 | |||
[交流]
线性基因组可视化神来之笔---karyoploteR
|
|
线性基因组可视化神来之笔---karyoploteR --光俊 随着人类基因组逐渐被破译，一张生命之图将被绘制，人们的生活也将发生巨大变化。今天我们给大家介绍一个包，对于性基因组可视化来说，可谓是神来之笔。 记得多年前，刚接触高通量测序，为了画一张覆盖度图，花了一周时间写代码。那时用的办法非常笨，不过总算能够完成操作。我们可以描述一下大概的过程：我们是从bam文件出发，先把bam文件转换为bed文件，有了bed文件以后，在对染色体划分区段，比如10kb为一段，然后计算每个区段的覆盖深度。然后在应用复杂的R代码，先把23挑染色体画出来，然后再把每个区段的深度plot上去。现在想想，听起来就觉得复杂。现在好了，这些复杂的工作，一个命令就搞定，这就是我们今天要介绍的karyoploteR包。 我们只要准备简单的输入文件，就能完成漂亮的基因图的绘制。比如我们做GEO和TCGA数据挖掘的时候，得到了差异表达的表格，我们只需要利用差异表达里面基因名字和logFC，就能完成一张漂亮的基因组图。 再比如，我们手头有CNV或者基因融合的结果，我们就可以利用最原始的数据，绘制基因组图。 更神奇的是，我们可以直接输入bam文件，就能得到漂亮的圈图。是不是很神奇，以前公司做项目，需要两三千的图形，现在只要几个命令就能搞定，想学的赶紧加入我们生信自学网。 下面我们看一段脚本，非常简单将bam直接可视化。 setwd("C:\\Users\\lexb4\\Desktop\\karyoploteR\\12.kpPlotBAMDensity"  library(karyoploteR) bamFile <- "input.bam" window.size <- 1e4 #compute the density with 10kb windows pdf(file="kpPlotBAMDensity.pdf",width=10,height=7) kp <- plotKaryotype(genome="hg38" kp <- kpPlotBAMDensity(kp, data = bamFile, window.size = window.size, r0=0, r1=1, ymax=600, col="darkorchid" kpAddBaseNumbers(kp, tick.dist=10000000, minor.tick.dist=1000000) dev.off() 接下来我们看看karyoploteR的效果图，如下： |
回复此楼» 猜你喜欢
欢迎发来filecode的Mz6后的代码验证其规律
已经有62人回复
-
时间戳又变了8-15
已经有23人回复
-
filecode=后面第一个是大写字母
已经有7人回复
-
有时候，自然基金真的不能太认真 （我的申报经验）
已经有10人回复
-
filecode
已经有12人回复
-
小木虫上这么多卖论文的，真有人买论文么？感觉没必要啊
已经有13人回复
-
哪位老哥知道今年的国自然具体哪一天放榜？
已经有7人回复
-
关于Filecode分析方法
已经有13人回复
-
各位道友，我要去昆明玩几天，回来见。
已经有8人回复
-
咱们一起用铁证分析2026国家社科基金中标与否
已经有22人回复
