ChIPseeker系列文章已经介绍了很多内容，包括注释的方方面面，也包括强大的可视化功能（《CS6: ChIPseeker的可视化方法（中秋节的视觉饕餮）》）。

今天要介绍一下数据挖掘，从大量已有的数据来产生新的hypothesis。正如我在ChIPseeker的文章里写的：

There are increasing evidences shown that combinations of TFs are important for regulating gene expression (Perez-Pinera et al., 2013; Zhu et al., 2008). However, systematically identification of TF interactions by ChIP-seq is still not available. Even if a specific TF binding is essential for a particular regulation was known, we do not have prior knowledge of all its co-factors. There are no systematic strategies available to identified un-known co-factors by ChIP- seq.

并没有方法可以大规模地预测未知的共同调控因子，而数据挖掘就是要给我们这种预测的能力。

我当年在写ChIPseeker的时候，我有纠结是写篇Bioinformatics的application note呢，还是写篇长文灌水NAR，毕竟NAR影响因子高一点，最后还是发了Bioinformatics，因为我没钱，囧，Bioinformatics不要版面费啊。然后限于篇幅，ChIPseeker有大量可视化的函数，我在文章中一张图都没放！！！如果当时决定发NAR的话，这个数据挖掘这一块我就会写多点。

Continue reading

meme无非是照片上加两行字，meme这个包本来是练手之作，用来模拟ggplot2的api。结果我后面又实现了字体阴影效果，这次我又带来了外部字体支持。

加载外部字体在R上面主要有extrafont和showtext两个实现，特别是showtext用起来还是非常简单的。我在meme包里写了font_import函数，内部就是用showtext，只是简单的包装了一下，这样小白用户可以直接下个字体就开撸。有阴影效果，有搞笑字体，就万事具备了，这俨然是个专业的meme工具了。

Continue reading

在《CS7：Genomic coordination的富集性分析（1）》说到了seq2pathway这个包，其实是两部曲，seq2gene->gene2pathway，无非是把测序片段用临近的基因注释，包括和TSS overlap的基因，宿主基因，上下游的基因等，然后拿这些基因跑ORA，做富集，仅此而已，这个包支持的物种极有限，《CS4：关于ChIPseq注释的几个问题》这一文中讲到ChIPseeker支持所有有基因组注释的物种，而《clusterProfiler for enrichment analysis》也支持所有物种（即使你自己跑的电子注释，也能支持），那么使用ChIPseeker来做基因注释，然后衔接clusterProfiler就可以支持所有物种的测序片段进行功能富集分析了。

《CS3: peak注释》本身就支持几种注释，另外我写了一个seq2gene的函数，套用seq2pathway的思路，把一个基因位置上所有关联的基因全部返回来，我们可以使用它去把基因位置信息转换成基因列表，然后用于富集分析，还是熟悉的味道，还是熟悉的配方🦄

Continue reading

Reassortment is an important strategy for influenza A viruses to introduce a HA subtype that is new to human populations, which creates the possibilities of pandemic.

A diagram showed above (Figure 2 of doi:10.1038/srep25549) is widely used to illustrate the reassortment events. While such diagrams are mostly manually draw and edit without software tool to automatically generate. Here, I implemented the hybrid_plot function for producing publication quality figure of reassortment events.

library(tibble)
library(ggplot2)


n <- 8

virus_info <- tibble(
    id = 1:7,
    x = c(rep(1990, 4), rep(2000, 2), 2009),
    y = c(1,2,3,5, 1.5, 3, 4),
    segment_color = list(
        rep('purple', n),
        rep('red', n),
        rep('darkgreen', n),
        rep('lightgreen', n),
        c('darkgreen', 'darkgreen', 'red', 'darkgreen', 'red', 'purple', 'red', 'purple'),
        c('darkgreen', 'darkgreen', 'red', 'darkgreen', 'darkgreen', 'purple', 'red', 'purple'),
        c('darkgreen', 'lightgreen', 'lightgreen', 'darkgreen', 'darkgreen', 'purple', 'red', 'purple'))
)

flow_info <- tibble(from = c(1,2,3,3,4,5,6),
                    to = c(5,5,5,6,7,6,7))
hybrid_plot(virus_info, flow_info)

Continue reading

取子集对于进化树可视化来说是非常常见的，我们要区分内部节点和外部节点，我们也可能想针对某些特点的节点进行注释。

ggplot2现在所有图层都不支持直接取子集，所以呢ggtree就自己定义了一些修改的图层，包括geom_text2, geom_label2, geom_point2 和 geom_segment2，这些图层和ggplot2的版本唯一差别就是支持取子集。这样对于我们做注释来说，就更方便了。

比如说我想给内部节点打点，可以用：

ggtree(tree) + geom_text2(aes(label=label, subset=!isTip), hjust=-.2) +
        geom_point2(aes(subset=!isTip), color="red", size=3)

Continue reading

我这个公众号不像大多数所谓的生信从入门到精通的各种其实只是搬运点入门教程的群众喜闻乐见的公众号。正如我在《为什么要开这个公众号》里说的，这是小众的，有个人色彩的各种原创文。我不可能像其它公众号一样招两小弟当客服，很多人在公众号后台向我扔了许多问题，由于个人精力有限，只能优先解答「知识星球」的问题，上次写的《同一数据多变量分组的boxplot?》，图虽然简单，却穿着好多件马甲，而我把它扒光了给你看🙈

上面这个图，你看着高大上吧，我都可以吐它一脸口水。每一个有灰色背景的图，在x轴上violin都够到边界了，其实所有的violin都够到了，这证明什么？每个violin之间其实不可比较！你能想像几个独立的数据，在统一的bin width情况下，画density curve，竟然最高点都一样高吗？显然可能性几乎为0。这个如果使用ggplot2的话，可以使用scale='width'强制拉成一样高，但我不推荐，正如我前面说的，不可比较了。默认参数scale='area'，积分面积一样，和density curve一样解析，另外的参数scale="count"，高度与计数同比例，和histogram一样解析，而scale='width'强制拉成一样高，如果没有在显眼处说明，误导性太强。

画这种图也可以手工拼，这样就简单了。在你需要的情况下，加个灰色背景嘛，最后拼图嘛。当然拼图不一定要在illustrator里拼，比如你用grid，先画好坐标轴，然后水平上定义几个一样大的viewport，每一个violin都画在相应的viewport里面，对于画图函数来说，viewport就是整个画布了（虽然只是画布里的一块区域），所以你要么画violin，要么在画之前先画个矩阵，一路画下来，代码可以直接生成这样的图，但这图每一个violin都是独立画的（当然也不是完全独立，每一个水平上的ylim是有统一的），就算代码一步生成，也跟illustrator里拼没两样。

这里我要教你用ggplot2自动生成，其实解决思路早已推送过，请看《facet_plot：加图层到特定分面，方法二》，也正如我在《什么！你的图上有一双看不见的手》里说的，你们以为我在教ggtree，其实同时在教ggplot2。

Continue reading

这个问题，其实答案就存在于《听说你还不会画热图》，我们先生成一个矩阵：

set.seed(2017-11-12)
d = data.frame(matrix(rnorm(100), ncol=10))
colnames(d) = paste0('t', 1:10)
rownames(d) = paste0('g', 1:10)

这个矩阵是rnorm生成的随机数，有正有负，我们再生成一个矩阵，只有正数，并且数值上比第一个矩阵要大：

d2 = abs(d) * 1.2

Continue reading