TADCompare：差异TAD分析

日常瞎掰

  上海这一波疫情来得真可谓猛烈！原本连“四五天”都不能封闭，属于全世界的上海都停摆了个把月时间。你说厉害不厉害！灾难来临时，趋利避害才是动物的天性！体现在人类身上就等价于人们常说的“生命健康至上”！作为一名安分守己的打工人，小编自然是听从政策的安排，从3月18号小区封闭以来就开始在家“坐月子”，刚开始觉得居家办公还挺好，时间长了不自在的感觉也就慢慢展现出来了。
  不止生活上打折扣，工作上也是诸多不便。庆幸的是13年买的第一台笔记本电脑没有罢工，这里真的要夸赞一下“华硕品质，坚如磐石”还真不是吹嘘的！不要小看这个“老年机”，在禁足的日子里，可是支撑起了两个人的工作内容，做到了从早上到晚上不间断地工作。不管如何，时间过得依然挺快，这一转眼已经“满月”了，期待再一转眼能够解封，回归以前正常的生活！

TADCompare

  最近一段时间一直在分析HiC-seq的数据，虽然对该技术有了一定的了解，但还是处于学习阶段。基因组折叠形成的拓扑结构域(topologically-associating domains, TADs)，在很多生物学过程如细胞分化、转录调控、基因组复制及DNA修复中有着重要的作用。研究TAD在不同处理中的变化也是HiC的用处之一。分析差异TAD的好用软件并不多，这里给大家分享一个——TADCompare，能不能得到大家的认可不知道，至少软件在考虑差异TAD方面考虑的还是挺全面的。
  该软件将TAD差异类型分为complex, split, merge, shifted, strength change，示意图如下：

  该软件定义差异TAD的过程还是挺有意思的。首先，基于Modified Spectral Clustering方法来计算交互矩阵的boundary score，然后用两个矩阵的score值计算得到一个gap_score值来评估是否为差异TAD。评估过程大致如下：

1. 如果满足gap_score > z_thresh(默认为2)，做为Differential，否则为Non-Differential，并且根据gap_score >0 来确定boundary score富集在matrix1或者matrix2，后续依据这个富集的情况再结合其他条件对Differential进行细致分类。
2. 接着将Differential中满足Bound_Dist(boundary的窗口滑动评估) <= 5划分为Shifted。
3. 接着将剩下的Differential的boundary按照前后分布为Non-Differential类型，若富集在matrix1中则为split，若富集在matrix2中则为merge。
4. 再将剩下的Differential的boundary中满足TAD_Score1 和 TAD_Score2 > 1.5条件的划为Strength Change类型。
5. 最后剩下的Differential的boundary即为Complex。PS：若提供两个矩阵的boundary文件来分析差异TAD，如果Boundary在提供的文件中出现且为Non-Differential则归类为Non-Overlap。
     上面的分类过程官网并没有提供，而是小编从软件的源码中整理所得，如有不对的地方希望大家莫见怪。下面我们来介绍一下软件的使用方法：

library(TADCompare)
library(dplyr)
library(SpectralTAD)

dim(rao_chr22_prim)
results = TADCompare(rao_chr22_prim, rao_chr22_rep, resolution = 50000)

#提供TAD
bed_coords1 = bind_rows(SpectralTAD::SpectralTAD(rao_chr22_prim, chr = "chr22", levels = 3))
bed_coords2 = bind_rows(SpectralTAD(rao_chr22_rep,  chr = "chr22", levels = 3))
Combined_Bed = list(bed_coords1, bed_coords2)
TD_Compare = TADCompare(rao_chr22_prim, rao_chr22_rep, resolution = 50000, pre_tads = Combined_Bed)

head(TD_Compare$TAD_Frame)

结果预览：

##   Boundary   Gap_Score TAD_Score1 TAD_Score2     Differential Enriched_In
## 1 16850000  0.23611002  7.8431724  7.5800558 Non-Differential    Matrix 1
## 2 17250000 -0.08875306  0.1103987  0.1409999 Non-Differential    Matrix 2
## 3 17600000 -0.06733403 -0.4153272 -0.3809658 Non-Differential    Matrix 2
## 4 18000000 -0.10338944  2.3616107  2.3473922 Non-Differential    Matrix 2
## 5 18350000 -0.04284963 -0.7767081 -0.7433986 Non-Differential    Matrix 2
## 6 19100000  0.01857312 -0.8269469 -0.8153997 Non-Differential    Matrix 1
##               Type
## 1      Non-Overlap
## 2      Non-Overlap
## 3      Non-Overlap
## 4 Non-Differential
## 5 Non-Differential
## 6      Non-Overlap

差异结果可视化：

p <- DiffPlot(tad_diff    = TD_Compare, 
              cont_mat1   = mtx1,
              cont_mat2   = mtx2,
              resolution  = res,
              start_coord = 8000000,
              end_coord   = 16000000,
              show_types  = TRUE, 
              point_size  = 5,
              max_height  = 5,
              rel_heights = c(1, 2),
              palette     = "RdYlBu")

plot(p)

结果如下：

  由上面的代码可知，使用过程并不复杂，这里不再详细解释。该软件支持多种matrix输入格式，如n×n、n×(n+3)、sparse 3-column、.hic files、.cool files、HiC-Pro files等矩阵格式。TADCompare的结果是一个包含三个内容的list，分别为TAD_Frame、Boundary_Scores、Count_Plot。前面两个为数据框，第一个是过滤后的结果，第二个包含全部结果，包含的列相同分别为Boundary、Gap_Score、TAD_Score1、TAD_Score2、Enriched_In、 Type。而Count_Plot是为展示每种类型数量的条形图的ggplot对象。

结束语

  上面的内容大概介绍了TADCompare主要功能，该软件还包含一些用于分析特定情况的函数如TimeCompare、ConsensusTADs，用法依然挺简单，具体适用于哪些情况，感兴趣的可以参考官网介绍：https://bioconductor.org/packages/devel/bioc/vignettes/TADCompare/inst/doc/TADCompare.html。