(PNAS 2019) scHiCluster （Part I：
关键词：single-cell, Hi-C，3D chromosomal structure，random walk
文章概要
本文的作者团队开发了scHiCluster 工具，使用linear convolution + RWR对single cell Hi-C数据进行Imputation。作者在文章中展示了imputation map的2个重要应用：（1）细胞聚类 （2）识别 TAD-like structure 。
原文：
Zhou J, Ma J, Chen Y, Cheng C, Bao B, Peng J, Sejnowski TJ, Dixon JR, Ecker JR. Robust single-cell Hi-C clustering by convolution- and random-walk-based imputation. Proc Natl Acad Sci U S A. 2019 Jul 9;116(28):14011-14018. doi: 10.1073/pnas.1901423116. Epub 2019 Jun 24. PMID: 31235599; PMCID: PMC6628819.
文章背景
自2013年起，多种单细胞Hi-C技术被发表（1-6），但是却缺乏用于单细胞Hi-C数据分析的计算工具，尤其是基于单细胞Hi-C数据进行细胞聚类。
在本文之前，基于单细胞Hi-C进行聚类的相关工作包括：
- 2018 Tan et al. （6）指出可以将在 bulk Hi-C 注释出的细胞特征用于 single-cell Hi-C数据，从而识别得到相应的细胞类型
- 2018 Liu et al.（7）提出 HiCRep + MDS 的embedding方法，将该方法用于2017, Nagano et al. 数据可以区分来自不同细胞周期的细胞
本文作者指出，单细胞Hi-C数据用于细胞分型的主要难点包括：
- 细胞内染色质三维结构具有高度的时空动态性，这导致了scHiC数据具有高易变性（Variability）
- 数据高度稀疏
- 数据深度（Coverage）往往成为驱动聚类的主要因素（Figure S1) ，但是不能保证导致不同细胞间Coverage变异度大的原因是实验因素，还是内在的生物学机制。
作者解决以上问题的思路是：
- 针对稀疏问题，作者使用Linear convolution + RWR （8）对数据首先进行imputation
- 针对高变异问题，作者完成Imputation后仅使用 top-ranked interactions
在2018年Liu et al. 的工作（HiCrep）中，作者在计算2个HiC map的相似性前首先使用 linear convolution 对raw contact matrice进行平滑处理； 而在2018年O. Ursu等人的工作（GenomeDISCO）（9）中，则使用首先使用 Random walk对数据进行平滑处理。 而本文的作者结合了2者的优点，同时将GenomeDISCO中使用的random walk 调整为 random walk with restart.
主要结果
Figure 1 | scHiCluster 原理
scHiCluster主要由4步构成）：
- 将Hi-C map中的每个单元格的值替换为该单元格与其周边单元格的加权平均数（linear convolution）
- 使用 random walk (with restart) 算法
- 仅保留top 20% interaction
- clustering & visualization
作者指出，convolution step 使三维互作信息在线性基因组上的近邻间传递，而之后的RWR 则使信息在网络/空间近邻间传递。
Figure 1
Figure 2 | 使用模拟数据集测试scHiCluster聚类效果
作者首先基于模拟数据测试算法效果，旨在测试算法在不同coverage以及不同resoluiton下的表现。
模拟数据生成
作者首先指出，如果直接对bulk Hi-C进行简单抽样（即downsample至相同contact数）得到的数据稀疏性低，变异度小，无法scHiC数据。因此，作者提出了一种特殊的抽样方法，通过对数据的稀疏性进行控制并人为地向数据中添加噪音，从而实现对单细胞HiC数据的模拟生成（Figure S2）。
作者将该方法应用于2个bulk Hi-C数据集：
| 年份 | 作者 | 细胞类型 | 参考文献 |
|---|---|---|---|
| 2014 | Rao et al. | GM12878, IMR90, HMEC, NHEK, K562, HUVEC, KBM7 | （10） |
| 2017 | Bonev et al. | mESC, NPC, CN | (11) |
作者测试了7种不同的coverage（500k, 250k, 100k, 50k, 25k, 10k, 5k）以及2种不同的resolution（1M, 200k）。在每种测试条件下，每种细胞类型模拟生成30个细胞的scHiC图谱。
结果评估标准
作者使用scHiCluster模拟数据进行无监督聚类，并使用校正的兰德系数（adjusted Rand index, ARI）评估分类的准确性。
评估结果
- 在2个数据集中，scHiCluster的表现均优于对照方法（直接使用PCA）（Figure S4）
- 当contacts < 25k时，scHiCluster性能下降； 5k contacts 时，无法聚类（Figure S5）
- 1M resolution 优于 200k（作者认为可能是因为分辨率更低时稀疏度也更低）
- window size的选择对结果影响不大（Figure S6）
Figure 3 | 在真实scHiC数据集测试scHiCluster聚类效果
用于测试的数据集
| 年份 | 作者 | 细胞类型 | # of contacts / cell | 参考文献 |
|---|---|---|---|---|
| 2017 | V. Ramani et al. | HeLa, HAP1, GM12878, K562 | 5.2k ~ 102.7k, median 97.3k | （2） |
| 2017 | Flyamer et al. | mouse zygotes and oocytes | 6.6k to 1.1M, median 97.3k | （5） |
| 2017 | Nagano et al. | mESC across different stages | （4） |
结果：
作者将scHiCluster与4种方法比较：1. PCA 2. HiCRep + MDS 3. eigenvector 4.Decay profile 。
结果无论是从可视化结果或使用ARI作为评估标准，scHiCluste均表现更优。（但是在Ramani数据集可视化结果中，scHiCluster仅Hela和其他细胞类型有明显区分，GM12878与K562有一定的区分度，但与HAP1混在一起）。
Figure 3 (A-B）
Restart probability与无监督聚类方法及参数的选择对结果的影响见Figure S8
与HiCRep+MDS方法相比，scHiCluster速度更快（Figure S9）。
对主成分的解读
PC1的weights matrix均平行于对角线，PC1更大的细胞中short-range contacts占比更高，作者认为这说明PC1可能捕捉到的是contact-distance curve信息，因此与细胞所处的细胞周期相关（Figure S10）。将scHiCluster应用于2017 Nagano数据集，来自不同周期的细胞可以被明显区分（Figure S11）。
Figure 4-5 | 鉴定 TAD-like structure 以注释细胞类型
作者首先举例说明，将来自同一细胞类型的不同细胞的 imputed contact matrices 整合到一起，可以观察到TAD样结构（TAD-like structures, TLSs） 。并且作者认为，不同细胞类型间的差异化TLS与差异化的基因表达相关，是注释细胞类型的关键（Figure S13-15）。
接下来，作者尝试使用 scHiCluster imputated contact matrices + TopDom （12）在全基因组范围识别TLS。
为了验证此方法的有效性，作者首先利用2017 Bonev et al.数据集，并选取了chr 19上一个特定位置（ESC和NPC的bulk Hi-C数据显示此位置上ESC和NPC细胞存在差异化的TAD）作者使用1 Mb resolution，基于不同coverage的数据进行TLS识别。结果显示：与raw contact matrices相比，imputed matrices 的可视化结果中TLS结构更为清晰，并且更有可能在单细胞水平上鉴定到差异化的TLS的边界。
Figure 4（A-B）
接下来作者利用Nagano et al.数据集证明了该方法在全基因组范围的有效性。
首先，作者将每个单细胞中检测到的TLS边界与bulk Hi-C的TAD boundary比较，平均46%的TLS boundary与TAD boundary重合，53%的TAD boundary在TLS boundary中能被检测到（Figure S16）。
接下来，作者统计并展示了1个10Mb区域内每个bin成为TLS边界的概率：
- 几乎每个bin都在至少一个单细胞中是TLS边界
- CTCF结合位点所在的bin有更高的概率是TLS边界
- TAD边界所在有更高的概率是TLS边界
Figure 5D
参考文献
[1] T. Nagano et al., Single-cell Hi-C reveals cell-to-cell variability in chromosome structure. Nature 502,59–64 (2013).
[2] V. Ramani et al., Massively multiplex single-cell Hi-C. Nat. Methods 14,263–266 (2017)
[3] T. J. Stevens et al., 3D structures of individual mammalian genomes studied by single- cell Hi-C. Nature 544,59–64 (2017).
[4] T. Nagano et al., (4)Cell-cycle dynamics of chromosomal organization at single-cell res- olution. Nature 547,61–67 (2017).
[5] I. M. Flyamer et al., Single-nucleus Hi-C reveals unique chromatin reorganization at oocyte-to-zygote transition. Nature 544, 110–114 (2017).
[6] L. Tan, D. Xing, C.-H. Chang, H. Li, X. S. Xie, Three-dimensional genome structures of single diploid human cells. Science 361, 924–928 (2018).
[7] J. Liu, D. Lin, G. G. Yardimci, W. S. Noble, Unsupervised embedding of single-cell Hi-C data. Bioinformatics 34,i96–i104 (2018).
[8] J.-Y. Pan, H.-J. Yang, C. Faloutsos, P. Duygulu, “Automatic multimedia cross-modal correlation discovery” in Proceedings of the Tenth ACM SIGKDD International Con- ference on Knowledge Discovery and Data Mining, KDD ’04 (ACM, New York, 2004), pp 653–658.
[9] O. Ursu et al., GenomeDISCO: A concordance score for chromosome conformation capture experiments using random walks on contact map graphs. Bioinformatics 34, 2701–2707 (2018).
[10] S. S. P. Rao et al., A 3D map of the human genome at kilobase resolution reveals principles of chromatin looping. Cell 159, 1665–1680 (2014).
[11] B. Bonev et al., Multiscale 3D genome rewiring during mouse neural development. Cell 171, 557–572.e24 (2017)
[12] H. Shin et al., TopDom: An efficient and deterministic method for identifying topo- logical domains in genomes. Nucleic Acids Res. 44, e70 (2016).
