染色質三維基因組結構在細胞功能和基因調控中起著關鍵作用。單細胞Hi-C技術可以在細胞水平上捕獲基因組結構信息，這為研究不同細胞類型之間基因組結構的變化提供了機會。最近，針對單細胞Hi-C數據分析出現了一些很好的計算分析方法。本文首先對可用的單細胞Hi-C數據分析方法進行綜述，包括單細胞Hi-C數據的預處理方法、基于單細胞Hi-C數據的多尺度結構識別方法、基于單細胞Hi-C數據集的類bulk Hi-C接觸矩陣生成方法、偽時間序列分析和細胞分類研究；然后闡述了單細胞Hi-C數據在細胞分化、結構變異的應用研究；最后展望了基于單細胞Hi-C數據的未來發展前景。
引用本文： 龔海燕, 麻付強, 張曉彤. 單細胞Hi-C數據分析及應用研究進展. 生物醫學工程學雜志, 2023, 40(5): 1033-1039. doi: 10.7507/1001-5515.202303046 復制
0 引言
真核生物的染色體在細胞核內如何折疊是一個古老的問題，至今仍未得到解答。早期的生化和基于顯微鏡技術的研究發現了染色質結構域和環，但這種組織特征的生物學意義尚不清楚。高通量染色體構象捕獲（high-throughput chromosome conformation capture，Hi-C）技術（例如BL-Hi-C[1]）提供了全基因組染色質交互作用，幫助研究者對間期核中染色體的組織有了新的認識，特別是，在廣泛的細胞類型中發現細胞類型不變、進化保守的拓撲關聯結構域[2]（topologically associated domain，TAD），為研究動物發育和人類疾病提供了新的分子框架。最近，單細胞技術允許研究者從單個細胞或單個細胞核中獲得Hi-C數據，從而識別單個細胞的染色質結構特征。單細胞Hi-C技術的觀察結果表明，TAD出現在細胞群，而不是作為靜態結構存在于單個細胞中。單細胞Hi-C技術為分析罕見細胞類型（如干細胞、腫瘤祖細胞、卵母細胞和全能性細胞）中的染色質折疊提供了機會，有助于更深入地了解發育和疾病的基本機制。
Hi-C技術的實驗步驟主要包括如下：① 通過甲醛交聯固定染色質，將由蛋白質介導的空間相鄰的染色質片段進行共價連接；② 使用限制性核酸內切酶切割基因組；③ 使用生物素標記切割的末端；④ 使用DNA連接酶連接切割的末端以產生DNA嵌合分子；⑤ 純化和破壞DNA嵌合分子，并篩選出具有生物素標簽的DNA片段；⑥ 對DNA文庫進行雙端測序。而單細胞Hi-C技術與傳統的基于群體細胞的Hi-C技術類似，主要的一個關鍵步驟是單個細胞的分離或者barcoding。典型的方法是將細胞或細胞核分離到單獨的反應混合物中，并分別執行后續步驟。例如，scHi-C技術[3]將傳統Hi-C協議的規模縮小到一個反應管一個細胞。
隨著單細胞Hi-C技術的發展，單細胞Hi-C數據經過數據質量判斷、消除單細胞Hi-C數據偏差等預處理操作，可進行以下研究：單細胞染色體三維結構建模、TAD識別、TAD之間關聯關系研究、基于單細胞Hi-C的染色體間交互作用研究、特異性染色質環研究、結構變異檢測、單細胞的偽時間序列、細胞之間的關系研究、區分細胞類型。近年來已有多位學者針對單細胞Hi-C數據的研究進展進行了綜述性分析，包括單細胞Hi-C協議的對比[4]、基于單細胞Hi-C數據的聚類方法[5-6]、細胞動態過程研究[7-8]、結構變異檢測工具[9]、單細胞中三維基因組組織現狀[10]、單細胞的整合方法綜述[11]等。本文從基于單細胞Hi-C數據計算方法的角度介紹了基于單細胞Hi-C數據分析的最新研究現狀，包括數據質量判斷、數據偏差消除、多尺度結構分析、偽時間序列分析以及細胞類別聚類分析。在此基礎上，本文分析了單細胞Hi-C數據的應用研究以及發展趨勢，作為現有單細胞Hi-C技術綜述成果在計算方面的更新和補充，以期為將三維基因組研究進一步拓展為四維基因組研究提供支持。
1 單細胞Hi-C預處理方法
如圖1a所示，單細胞Hi-C技術的一個重要步驟就是對單個細胞進行Hi-C測序，或者通過barcoding的方式得到單細胞Hi-C測序讀對數據。由于單細胞檢測技術固有的高水平噪聲阻礙了進一步研究，因此，在進行單細胞Hi-C數據分析之前，需對單細胞Hi-C數據進行質量評估。如附件1所示，2020年，Horton等[12]提出工具GiniQC，通過計算順勢讀對百分比、基尼指數來衡量數據質量，GiniQC值越高，代表數據質量越高。Wolff等[13]在Galaxy HiCExplorer中集成了hicQuickQC，只需要使用一部分數量的讀對就可以計算單細胞Hi-C數據質量。特別地，使用hicQuickQC為前100萬次讀取創建質量報告的計算時間小于3 min。
圖1
基于單細胞Hi-C的數據分析方法及應用
a. 單細胞Hi-C數據獲取示意圖；b. 基于單個單細胞Hi-C數據的分析；c. 基于單細胞Hi-C數據集合的分析與應用
Figure1. Data analysis method and application based on single-cell Hi-Ca. schematic diagram of single-cell Hi-C data acquisition; b. analysis based on single single-cell Hi-C data; c. analysis and application based on bulk single-cell Hi-C datasets
在進行下游分析前，需對得到的單細胞Hi-C測序讀對數據進行以下預處理過程：讀對預處理（可選）、讀對映射（可選）、限制性片段分配（可選）、接觸值過濾、重復數據刪除，最后通過分位點生成單細胞Hi-C接觸矩陣。單細胞Hi-C接觸矩陣是一個對稱的位點對位點矩陣，每個條目表示潛在空間接近的位點對之間的相互作用頻率值。由于單細胞Hi-C數據容易出現技術噪聲和偏差，阻礙下游分析，因此，需對單細胞Hi-C接觸矩陣進行去噪和歸一化后進行下游分析。2022年，Zheng等[14]提出一種快速縮放歸一化方法BandNorm，以解釋scHi-C特定的偏差。BandNorm方法首先去除細胞內的基因組距離偏差，然后在細胞之間進行測序深度歸一化，再為細胞間接觸矩陣加上帶依賴接觸衰減估計值，最終得到歸一化后接觸矩陣。在基準測試實驗中，BandNorm在細胞類型分離、相互作用位點識別和細胞類型關系恢復方面表現領先。
經過數據預處理得到的單細胞Hi-C接觸矩陣主要以scool格式[15]進行存儲，不僅高效、用戶友好、節省存儲，還可保證穩定的API支持。后續單細胞Hi-C數據分析主要分為兩種，一種為針對單個單細胞Hi-C接觸矩陣的數據分析，如圖1b所示，主要是基于單細胞的多尺度結構識別方法（包括單細胞三維結構建模、單細胞A/B區室識別、單細胞TAD識別、單細胞Loop識別）。另一種分析如圖1c所示，針對包含多個細胞的單細胞Hi-C數據集的數據分析，主要包括偽時間序列分析、通過聚類或者嵌入的方法進行細胞聚類分析或通過插補技術得到一個類bulk Hi-C接觸矩陣后進行平均TAD或Loop識別。在下文中，將主要針對兩種分析類型和基于單細胞Hi-C數據的應用研究進行闡述。
2 基于單細胞Hi-C數據的多尺度結構識別方法
在基于群體細胞得到的Hi-C數據分析中，可從染色質三維結構[16]、A/B區室[17]、TAD[18]和Loop[19]四個尺度進行識別。同樣的，如附件2所示，目前該領域也研發出針對單個單細胞Hi-C接觸矩陣的算法，進行單細胞三維結構建模、單細胞A/B區室識別、單細胞TAD識別以及單細胞Loop識別。
從染色質三維結構建模的尺度來看，2021年Kos等[20]基于接觸位點的彈簧松弛特性提出聚合物模型DPDchrom。Meng等[21]基于貝葉斯理論框架開發了單細胞染色體構象計算工具Si-C，實現了10 kb分辨率下的單細胞三維結構重構。Messelink等[22]提出一種完全數據驅動的、基于最大熵的計算方法MaxEnt，從模型生物新月形的單細胞Hi-C數據中重構了三維染色體構象。Zha等[23]利用Lennard-Jones勢、Metropolis-Hastings算法對不同位置的DNA鏈珠進行測試，采用模擬退火算法對損失函數進行優化，重建了基于500 kb和50 kb分辨率單細胞Hi-C數據的染色體三維結構。但這些方法仍然是數據驅動的，并沒有解釋染色質結構形成的實際機制[24]。
有兩種方法可以從單細胞Hi-C接觸矩陣中識別A/B區室、TAD和Loop。一種是直接通過特征提取識別；一種是通過聚合的方式，得到一個偽bulk Hi-C接觸矩陣后，再用已有的bulk Hi-C方法進行識別。本節主要針對特征提取的方式進行闡述。在A/B區室識別方法上，Polovnikov等[25]基于聚合物非回溯算子的光譜方法，根據最大相關熵原則解析區室結構域，發現在單細胞Hi-C中識別的染色質區室具有生物學意義。在TAD識別方法上，Li等[26]提出deTOKI方法，利用非負矩陣分解尋找將基因組隔離成具有最小聚類塊的區域。2023年，Ye等[27]提出了一種分層染色質域結構識別算法（HiCS），將分層染色質結構域的識別問題轉化為在不同基因組尺度上尋找絕緣強度的峰值。在Loop識別方法上，2022年Zhang等[28]提供了DeepLoop工具，通過深度學習實現對Hi-C接觸矩陣的信號增強，以從低測序深度的Hi-C數據（包括單細胞Hi-C數據）中識別Loop。
3 基于單細胞Hi-C數據集的分析方法
3.1 類bulk Hi-C接觸矩陣生成方法
由于單細胞Hi-C接觸矩陣稀疏，無法直接使用現有的bulk Hi-C的計算方法進行計算，因此，需先對單細胞Hi-C數據進行插值處理。如附件3所示，2019年Zhou等[29]提出scHiCluster方法，將染色質交互作用看作一個網絡，首先采用平均的方法對Hi-C接觸值進行卷積，然后使用隨機游走對接觸值進行平滑，從而得到插值后的Hi-C接觸矩陣，降低了單細胞Hi-C接觸矩陣的稀疏性。Liu等[30]提出scHiCEmbed方法，從單細胞Hi-C數據中學習節點嵌入。scHiCEmbed將一個原始的或經schicluster估算的單細胞Hi-C接觸矩陣作為圖數據的鄰接矩陣，通過學習過程來獲得單細胞Hi-C矩陣的潛在表示，這些潛在表示即作為插值后矩陣。Xie等[31]提出了一種貝葉斯層次模型HiCImpute，通過考慮單細胞Hi-C二維數據結構的空間依賴性，同時也從類似的單細胞和批量Hi-C數據中借用信息，進行數據插值，從而提高數據質量。2023年，Liu等[32]提出scDEC-Hi-C方法，借助譯碼器作為插值后的單細胞Hi-C接觸矩陣。結合近兩年的插補方法可以看出，隨著圖神經網絡的發展，加上Hi-C接觸矩陣具備的圖特性，利用基于圖的編碼器-譯碼器模型可有效對單細胞Hi-C接觸矩陣插值。
3.2 基于單細胞Hi-C的偽時間序列分析
單細胞Hi-C技術為高分辨率下闡明染色體動力學提供了前所未有的機會。如何利用單細胞Hi-C接觸矩陣構建偽時間路徑對細胞沿發育軌跡進行排序是一個具有挑戰性的課題，因為該技術產生的這些矩陣具有高維和稀疏性，且受到噪聲和偏差的影響，并且它們所對應的軌跡拓撲可能是多樣化的。2019年，Ye等[33]開發了工具CIRCLET，通過考慮染色體結構的多尺度特征，在不指定起始細胞的情況下，利用1 171個單細胞Hi-C數據解析了單細胞的細胞周期階段，將軌跡劃分為12個階段，有助于揭示與動態亞結構相關的重要調控基因，為在單細胞分辨率下發現調控區域甚至癌癥標志物提供了新的框架。CIRCLET算法入選“2019年度中國生物信息學十大算法”。2021年，Bonora等[34]通過整合高通量單細胞RNA-seq、ATAC-seq和Hi-C獲得這三種模式的等位基因特異性數據，首次追蹤了與小鼠胚胎干細胞分化和X失活相關的變化，利用單細胞Hi-C技術獲得的等位基因特異性接觸衰減圖譜發現，失活的X染色體在X失活的分化細胞中具有獨特的圖譜。2022年，Lyu等[35]提出一種無監督的基于圖的工具scHiCPTR，用于從單細胞Hi-C接觸矩陣推斷具有生物學意義的偽時間軌跡。scHiCPTR提供了一個包括嵌入、圖構造、對偶圖細化、偽時間計算和結果可視化的工作流程。scHiCPTR能夠處理具有多種拓撲結構的發育軌跡，包括線性、分叉和圓形，并且與針對單細胞RNA-seq數據開發的方法具有競爭力。雖然已有上述工具可以進行偽時間序列分析，但是細胞發育過程中的軌跡拓撲可能是多樣化的，包括閉合環形、線性、單分叉、多分叉、樹形、連接圖、斷開圖等多樣性結構，不同的方法適用于不同的軌跡拓撲結構，而且受軌跡拓撲結構影響較大。如何在考慮不同拓撲結構的情況下，準確、高效地進行偽時間路徑構建和排序，是一個具有挑戰性的問題。
3.3 基于單細胞Hi-C數據集的細胞類別聚類分析
為了研究不同細胞類型之間染色體結構的變化，需要利用稀疏、異構的單細胞Hi-C數據計算方法。然而，很少有方法能夠準確有效地將這些數據按照細胞類型進行聚類。如附件4所示，2019年Zhou等[29]提出了一種針對單細胞Hi-C接觸矩陣的聚類算法scHiCluster。scHiCluster基于線性卷積和隨機游走進行插值，以模擬的和真實的單細胞Hi-C數據為基準，在低覆蓋率數據集上顯著提高了聚類精度。scHiCluster方法可在單個細胞內識別出TAD，但它缺乏捕獲可區分細胞類型的染色質交互作用的能力。2020年，Kim等[36]采用一種針對文檔集合的主題建模方法潛在狄利克雷分配（latent Dirichlet allocation，LDA），假設每個細胞對應一個文檔，而每個單細胞Hi-C交互作用對應一個單詞，揭示了染色質結構特征和特定細胞類型特征之間的關系。
單細胞Hi-C接觸矩陣非常高維、稀疏。數千個單細胞Hi-C接觸矩陣被聚合、展平為一個矩陣。根據分辨率的不同，該矩陣可能有幾百萬甚至幾十億個特征，因此，計算可能需要占用大量內存。2021年，Wolff等[37]提出了一種基于局部敏感哈希函數MINHASH 實現近似最近鄰的單細胞Hi-C聚類方法scHiCExplorer，通過計算K-最近鄰圖來減少相對于細胞數目的特征，從而降低維度和計算資源。隨著超圖技術的發展，Higashi算法[38]將單細胞Hi-C數據轉換為一張超圖，其中每個細胞和每個基因組位點分別表示為細胞節點和基因組位點節點，而單細胞接觸矩陣中的每個非零值被建模為連接相應細胞節點和該特定染色質交互作用的兩個基因組位點的超邊。Higashi利用基于超圖表示學習的算法，允許跨細胞協調信息，從而利用細胞之間的潛在相關性進行嵌入。2022年，Zhang等[39]提出一種基于張量分解和部分隨機游走的超快和可解釋的方法Fast-Higashi。Fast-Higashi能夠從稀疏的單細胞Hi-C數據中識別細胞類型，可以直接識別定義不同細胞類型的三維基因組特征，并幫助闡明基因組結構和功能之間的細胞類型特異性連接。Wu等[40]提出了一種基于單細胞Hi-C數據的高精度細胞分類算法scHiCStackL。scHiCStackL方法首先改進了單細胞Hi-C數據預處理方法，使生成的嵌入可更好地表示細胞，然后通過構造一個兩層堆疊集成模型進行細胞分類。2023年，Liu等[32]提出一種基于深度生成神經網絡的單細胞Hi-C數據分析框架scDEC-Hi-C，可揭示不同細胞類型的染色質結構差異。Zhou等[41]提出了一個高性能的單細胞Hi-C聚類框架 scHiCSC。作者首先提出了一種新的基于接觸數權重的平滑方法，以生成更精確的細胞特征嵌入；然后提出了一種新的特征融合方法，通過融合細胞內染色體結構信息和細胞間距離信息，進一步補充細胞的特征信息。實驗結果表明，scHiCSC對不同規模的數據集具有較強的泛化能力，性能優于現有的單細胞Hi-C聚類框架。雖然已經有針對單細胞Hi-C數據的線性或者非線性嵌入技術，但不同于基于單細胞RNA-seq數據的嵌入方法，目前尚未出現基于單細胞Hi-C數據的流形學習技術。單細胞Hi-C可以以矩陣、圖、向量的不同形式進行表達，這允許研究者們基于這些數據集進行更多的數據嵌入研究。
4 基于單細胞Hi-C數據的應用研究
對于器官/組織發育而言，構建人體細胞圖譜是了解器官功能的基礎[42]。單細胞組學在腦科學、腫瘤異質性、腫瘤免疫學、腫瘤耐藥和藥物開發等方面有一定的應用[43-44]。胚胎干細胞可以通過將自身逐漸暴露于特定因子，使用譜系特異性基因表達程序來產生功能細胞類型。Khateb等[45]使用單細胞Hi-C、RNA-seq、ATAC-seq等進行了胚胎干細胞多能性的染色質構像、單/體細胞基因表達、組蛋白修飾和可及性分析，獲得了前體細胞前中胚層的fate，并進一步進行了肌源性和神經源性分化，揭示了基因組和轉錄特征，并確定了指導初始Pax7表達和肌源性和神經源性過程激活的調控區域。Wang等[46]開發了EagleC工具，在慢性髓系白血病細胞系HAP1和K562上可預測染色體間移位的結構變異。也就是說，基于單細胞Hi-C數據，可進行基因結構變異的檢測，從而進行基因疾病的治療。
免疫細胞在分化和免疫應答過程中動態地重新編程其基因表達譜。然而，潛在的機制仍然難以捉摸。Lin等[47]開發了一種單細胞Hi-C方法，系統地描繪了巨噬細胞在這些過程中的三維基因組和動態表觀遺傳圖譜，并提出了“無序度”來測量拓撲關聯結構域內的基因組組織模式，這與單個細胞中的染色質表觀遺傳狀態、基因表達和染色質結構變異性相關。作者還發現NF-kappa B在結核分枝桿菌感染時啟動了系統性染色質構象重組。因此，單細胞Hi-C技術可以為研究細胞發育、細胞分化周期提供結構信息。
5 總結與展望
單細胞分辨率下分析三維基因組，為研究高階染色質結構和功能的異質性和動態性提供了獨特的機會。然而，單細胞Hi-C數據集的稀疏性、噪聲性和高維性等問題也帶來了新的分析挑戰。正如我們在這篇綜述中所討論的，計算方法已經專門用于單細胞Hi-C數據的分析，包括數據處理、降維、接觸圖插補、三維基因組結構建模和其他下游分析（例如，在不同尺度下調用三維基因組特征）。迄今為止，大多數單細胞Hi-C數據的計算方法都集中在降維和多尺度結構建模上，而其他類型的分析和解釋的重要方向需要更多的工作。
除了專門設計用于分析單細胞Hi-C數據的計算方法外，開發能夠將單細胞Hi-C數據與其他單細胞數據，例如單細胞RNA-seq數據、單細胞ATAC-seq數據等進行聯合分析的框架也至關重要。例如，聯合分析不同單細胞基因組數據的方法有望極大地促進對三維基因組結構、表觀基因組特征和轉錄活動的系統級理解。然而，單細胞Hi-C數據的極端稀疏性以及缺乏將單細胞Hi-C數據與其他單細胞分析相結合的系統評估，使得可靠的多模式分析具有相當的挑戰性。一個能同時分析單細胞三維基因組和其他表觀基因組以及功能特征的計算工具，將提供單細胞基因組結構和功能之間的直接關聯。總之，作為三維基因組和表觀基因組學研究中令人振奮的新方向，單細胞Hi-C和其他單細胞技術有望為核基因組結構和功能提供獨特和新穎的見解，這也將為開發新的計算方法創造機會。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：龔海燕查閱、整理相關文獻資料，麻付強起草論文，張曉彤負責論文整體框架和校對。
本文附件見本刊網站的電子版本（biomedeng.cn）。
0 引言
真核生物的染色體在細胞核內如何折疊是一個古老的問題，至今仍未得到解答。早期的生化和基于顯微鏡技術的研究發現了染色質結構域和環，但這種組織特征的生物學意義尚不清楚。高通量染色體構象捕獲（high-throughput chromosome conformation capture，Hi-C）技術（例如BL-Hi-C[1]）提供了全基因組染色質交互作用，幫助研究者對間期核中染色體的組織有了新的認識，特別是，在廣泛的細胞類型中發現細胞類型不變、進化保守的拓撲關聯結構域[2]（topologically associated domain，TAD），為研究動物發育和人類疾病提供了新的分子框架。最近，單細胞技術允許研究者從單個細胞或單個細胞核中獲得Hi-C數據，從而識別單個細胞的染色質結構特征。單細胞Hi-C技術的觀察結果表明，TAD出現在細胞群，而不是作為靜態結構存在于單個細胞中。單細胞Hi-C技術為分析罕見細胞類型（如干細胞、腫瘤祖細胞、卵母細胞和全能性細胞）中的染色質折疊提供了機會，有助于更深入地了解發育和疾病的基本機制。
Hi-C技術的實驗步驟主要包括如下：① 通過甲醛交聯固定染色質，將由蛋白質介導的空間相鄰的染色質片段進行共價連接；② 使用限制性核酸內切酶切割基因組；③ 使用生物素標記切割的末端；④ 使用DNA連接酶連接切割的末端以產生DNA嵌合分子；⑤ 純化和破壞DNA嵌合分子，并篩選出具有生物素標簽的DNA片段；⑥ 對DNA文庫進行雙端測序。而單細胞Hi-C技術與傳統的基于群體細胞的Hi-C技術類似，主要的一個關鍵步驟是單個細胞的分離或者barcoding。典型的方法是將細胞或細胞核分離到單獨的反應混合物中，并分別執行后續步驟。例如，scHi-C技術[3]將傳統Hi-C協議的規模縮小到一個反應管一個細胞。
隨著單細胞Hi-C技術的發展，單細胞Hi-C數據經過數據質量判斷、消除單細胞Hi-C數據偏差等預處理操作，可進行以下研究：單細胞染色體三維結構建模、TAD識別、TAD之間關聯關系研究、基于單細胞Hi-C的染色體間交互作用研究、特異性染色質環研究、結構變異檢測、單細胞的偽時間序列、細胞之間的關系研究、區分細胞類型。近年來已有多位學者針對單細胞Hi-C數據的研究進展進行了綜述性分析，包括單細胞Hi-C協議的對比[4]、基于單細胞Hi-C數據的聚類方法[5-6]、細胞動態過程研究[7-8]、結構變異檢測工具[9]、單細胞中三維基因組組織現狀[10]、單細胞的整合方法綜述[11]等。本文從基于單細胞Hi-C數據計算方法的角度介紹了基于單細胞Hi-C數據分析的最新研究現狀，包括數據質量判斷、數據偏差消除、多尺度結構分析、偽時間序列分析以及細胞類別聚類分析。在此基礎上，本文分析了單細胞Hi-C數據的應用研究以及發展趨勢，作為現有單細胞Hi-C技術綜述成果在計算方面的更新和補充，以期為將三維基因組研究進一步拓展為四維基因組研究提供支持。
1 單細胞Hi-C預處理方法
如圖1a所示，單細胞Hi-C技術的一個重要步驟就是對單個細胞進行Hi-C測序，或者通過barcoding的方式得到單細胞Hi-C測序讀對數據。由于單細胞檢測技術固有的高水平噪聲阻礙了進一步研究，因此，在進行單細胞Hi-C數據分析之前，需對單細胞Hi-C數據進行質量評估。如附件1所示，2020年，Horton等[12]提出工具GiniQC，通過計算順勢讀對百分比、基尼指數來衡量數據質量，GiniQC值越高，代表數據質量越高。Wolff等[13]在Galaxy HiCExplorer中集成了hicQuickQC，只需要使用一部分數量的讀對就可以計算單細胞Hi-C數據質量。特別地，使用hicQuickQC為前100萬次讀取創建質量報告的計算時間小于3 min。
圖1
基于單細胞Hi-C的數據分析方法及應用
a. 單細胞Hi-C數據獲取示意圖；b. 基于單個單細胞Hi-C數據的分析；c. 基于單細胞Hi-C數據集合的分析與應用
Figure1. Data analysis method and application based on single-cell Hi-Ca. schematic diagram of single-cell Hi-C data acquisition; b. analysis based on single single-cell Hi-C data; c. analysis and application based on bulk single-cell Hi-C datasets
在進行下游分析前，需對得到的單細胞Hi-C測序讀對數據進行以下預處理過程：讀對預處理（可選）、讀對映射（可選）、限制性片段分配（可選）、接觸值過濾、重復數據刪除，最后通過分位點生成單細胞Hi-C接觸矩陣。單細胞Hi-C接觸矩陣是一個對稱的位點對位點矩陣，每個條目表示潛在空間接近的位點對之間的相互作用頻率值。由于單細胞Hi-C數據容易出現技術噪聲和偏差，阻礙下游分析，因此，需對單細胞Hi-C接觸矩陣進行去噪和歸一化后進行下游分析。2022年，Zheng等[14]提出一種快速縮放歸一化方法BandNorm，以解釋scHi-C特定的偏差。BandNorm方法首先去除細胞內的基因組距離偏差，然后在細胞之間進行測序深度歸一化，再為細胞間接觸矩陣加上帶依賴接觸衰減估計值，最終得到歸一化后接觸矩陣。在基準測試實驗中，BandNorm在細胞類型分離、相互作用位點識別和細胞類型關系恢復方面表現領先。
經過數據預處理得到的單細胞Hi-C接觸矩陣主要以scool格式[15]進行存儲，不僅高效、用戶友好、節省存儲，還可保證穩定的API支持。后續單細胞Hi-C數據分析主要分為兩種，一種為針對單個單細胞Hi-C接觸矩陣的數據分析，如圖1b所示，主要是基于單細胞的多尺度結構識別方法（包括單細胞三維結構建模、單細胞A/B區室識別、單細胞TAD識別、單細胞Loop識別）。另一種分析如圖1c所示，針對包含多個細胞的單細胞Hi-C數據集的數據分析，主要包括偽時間序列分析、通過聚類或者嵌入的方法進行細胞聚類分析或通過插補技術得到一個類bulk Hi-C接觸矩陣后進行平均TAD或Loop識別。在下文中，將主要針對兩種分析類型和基于單細胞Hi-C數據的應用研究進行闡述。
2 基于單細胞Hi-C數據的多尺度結構識別方法
在基于群體細胞得到的Hi-C數據分析中，可從染色質三維結構[16]、A/B區室[17]、TAD[18]和Loop[19]四個尺度進行識別。同樣的，如附件2所示，目前該領域也研發出針對單個單細胞Hi-C接觸矩陣的算法，進行單細胞三維結構建模、單細胞A/B區室識別、單細胞TAD識別以及單細胞Loop識別。
從染色質三維結構建模的尺度來看，2021年Kos等[20]基于接觸位點的彈簧松弛特性提出聚合物模型DPDchrom。Meng等[21]基于貝葉斯理論框架開發了單細胞染色體構象計算工具Si-C，實現了10 kb分辨率下的單細胞三維結構重構。Messelink等[22]提出一種完全數據驅動的、基于最大熵的計算方法MaxEnt，從模型生物新月形的單細胞Hi-C數據中重構了三維染色體構象。Zha等[23]利用Lennard-Jones勢、Metropolis-Hastings算法對不同位置的DNA鏈珠進行測試，采用模擬退火算法對損失函數進行優化，重建了基于500 kb和50 kb分辨率單細胞Hi-C數據的染色體三維結構。但這些方法仍然是數據驅動的，并沒有解釋染色質結構形成的實際機制[24]。
有兩種方法可以從單細胞Hi-C接觸矩陣中識別A/B區室、TAD和Loop。一種是直接通過特征提取識別；一種是通過聚合的方式，得到一個偽bulk Hi-C接觸矩陣后，再用已有的bulk Hi-C方法進行識別。本節主要針對特征提取的方式進行闡述。在A/B區室識別方法上，Polovnikov等[25]基于聚合物非回溯算子的光譜方法，根據最大相關熵原則解析區室結構域，發現在單細胞Hi-C中識別的染色質區室具有生物學意義。在TAD識別方法上，Li等[26]提出deTOKI方法，利用非負矩陣分解尋找將基因組隔離成具有最小聚類塊的區域。2023年，Ye等[27]提出了一種分層染色質域結構識別算法（HiCS），將分層染色質結構域的識別問題轉化為在不同基因組尺度上尋找絕緣強度的峰值。在Loop識別方法上，2022年Zhang等[28]提供了DeepLoop工具，通過深度學習實現對Hi-C接觸矩陣的信號增強，以從低測序深度的Hi-C數據（包括單細胞Hi-C數據）中識別Loop。
3 基于單細胞Hi-C數據集的分析方法
3.1 類bulk Hi-C接觸矩陣生成方法
由于單細胞Hi-C接觸矩陣稀疏，無法直接使用現有的bulk Hi-C的計算方法進行計算，因此，需先對單細胞Hi-C數據進行插值處理。如附件3所示，2019年Zhou等[29]提出scHiCluster方法，將染色質交互作用看作一個網絡，首先采用平均的方法對Hi-C接觸值進行卷積，然后使用隨機游走對接觸值進行平滑，從而得到插值后的Hi-C接觸矩陣，降低了單細胞Hi-C接觸矩陣的稀疏性。Liu等[30]提出scHiCEmbed方法，從單細胞Hi-C數據中學習節點嵌入。scHiCEmbed將一個原始的或經schicluster估算的單細胞Hi-C接觸矩陣作為圖數據的鄰接矩陣，通過學習過程來獲得單細胞Hi-C矩陣的潛在表示，這些潛在表示即作為插值后矩陣。Xie等[31]提出了一種貝葉斯層次模型HiCImpute，通過考慮單細胞Hi-C二維數據結構的空間依賴性，同時也從類似的單細胞和批量Hi-C數據中借用信息，進行數據插值，從而提高數據質量。2023年，Liu等[32]提出scDEC-Hi-C方法，借助譯碼器作為插值后的單細胞Hi-C接觸矩陣。結合近兩年的插補方法可以看出，隨著圖神經網絡的發展，加上Hi-C接觸矩陣具備的圖特性，利用基于圖的編碼器-譯碼器模型可有效對單細胞Hi-C接觸矩陣插值。
3.2 基于單細胞Hi-C的偽時間序列分析
單細胞Hi-C技術為高分辨率下闡明染色體動力學提供了前所未有的機會。如何利用單細胞Hi-C接觸矩陣構建偽時間路徑對細胞沿發育軌跡進行排序是一個具有挑戰性的課題，因為該技術產生的這些矩陣具有高維和稀疏性，且受到噪聲和偏差的影響，并且它們所對應的軌跡拓撲可能是多樣化的。2019年，Ye等[33]開發了工具CIRCLET，通過考慮染色體結構的多尺度特征，在不指定起始細胞的情況下，利用1 171個單細胞Hi-C數據解析了單細胞的細胞周期階段，將軌跡劃分為12個階段，有助于揭示與動態亞結構相關的重要調控基因，為在單細胞分辨率下發現調控區域甚至癌癥標志物提供了新的框架。CIRCLET算法入選“2019年度中國生物信息學十大算法”。2021年，Bonora等[34]通過整合高通量單細胞RNA-seq、ATAC-seq和Hi-C獲得這三種模式的等位基因特異性數據，首次追蹤了與小鼠胚胎干細胞分化和X失活相關的變化，利用單細胞Hi-C技術獲得的等位基因特異性接觸衰減圖譜發現，失活的X染色體在X失活的分化細胞中具有獨特的圖譜。2022年，Lyu等[35]提出一種無監督的基于圖的工具scHiCPTR，用于從單細胞Hi-C接觸矩陣推斷具有生物學意義的偽時間軌跡。scHiCPTR提供了一個包括嵌入、圖構造、對偶圖細化、偽時間計算和結果可視化的工作流程。scHiCPTR能夠處理具有多種拓撲結構的發育軌跡，包括線性、分叉和圓形，并且與針對單細胞RNA-seq數據開發的方法具有競爭力。雖然已有上述工具可以進行偽時間序列分析，但是細胞發育過程中的軌跡拓撲可能是多樣化的，包括閉合環形、線性、單分叉、多分叉、樹形、連接圖、斷開圖等多樣性結構，不同的方法適用于不同的軌跡拓撲結構，而且受軌跡拓撲結構影響較大。如何在考慮不同拓撲結構的情況下，準確、高效地進行偽時間路徑構建和排序，是一個具有挑戰性的問題。
3.3 基于單細胞Hi-C數據集的細胞類別聚類分析
為了研究不同細胞類型之間染色體結構的變化，需要利用稀疏、異構的單細胞Hi-C數據計算方法。然而，很少有方法能夠準確有效地將這些數據按照細胞類型進行聚類。如附件4所示，2019年Zhou等[29]提出了一種針對單細胞Hi-C接觸矩陣的聚類算法scHiCluster。scHiCluster基于線性卷積和隨機游走進行插值，以模擬的和真實的單細胞Hi-C數據為基準，在低覆蓋率數據集上顯著提高了聚類精度。scHiCluster方法可在單個細胞內識別出TAD，但它缺乏捕獲可區分細胞類型的染色質交互作用的能力。2020年，Kim等[36]采用一種針對文檔集合的主題建模方法潛在狄利克雷分配（latent Dirichlet allocation，LDA），假設每個細胞對應一個文檔，而每個單細胞Hi-C交互作用對應一個單詞，揭示了染色質結構特征和特定細胞類型特征之間的關系。
單細胞Hi-C接觸矩陣非常高維、稀疏。數千個單細胞Hi-C接觸矩陣被聚合、展平為一個矩陣。根據分辨率的不同，該矩陣可能有幾百萬甚至幾十億個特征，因此，計算可能需要占用大量內存。2021年，Wolff等[37]提出了一種基于局部敏感哈希函數MINHASH 實現近似最近鄰的單細胞Hi-C聚類方法scHiCExplorer，通過計算K-最近鄰圖來減少相對于細胞數目的特征，從而降低維度和計算資源。隨著超圖技術的發展，Higashi算法[38]將單細胞Hi-C數據轉換為一張超圖，其中每個細胞和每個基因組位點分別表示為細胞節點和基因組位點節點，而單細胞接觸矩陣中的每個非零值被建模為連接相應細胞節點和該特定染色質交互作用的兩個基因組位點的超邊。Higashi利用基于超圖表示學習的算法，允許跨細胞協調信息，從而利用細胞之間的潛在相關性進行嵌入。2022年，Zhang等[39]提出一種基于張量分解和部分隨機游走的超快和可解釋的方法Fast-Higashi。Fast-Higashi能夠從稀疏的單細胞Hi-C數據中識別細胞類型，可以直接識別定義不同細胞類型的三維基因組特征，并幫助闡明基因組結構和功能之間的細胞類型特異性連接。Wu等[40]提出了一種基于單細胞Hi-C數據的高精度細胞分類算法scHiCStackL。scHiCStackL方法首先改進了單細胞Hi-C數據預處理方法，使生成的嵌入可更好地表示細胞，然后通過構造一個兩層堆疊集成模型進行細胞分類。2023年，Liu等[32]提出一種基于深度生成神經網絡的單細胞Hi-C數據分析框架scDEC-Hi-C，可揭示不同細胞類型的染色質結構差異。Zhou等[41]提出了一個高性能的單細胞Hi-C聚類框架 scHiCSC。作者首先提出了一種新的基于接觸數權重的平滑方法，以生成更精確的細胞特征嵌入；然后提出了一種新的特征融合方法，通過融合細胞內染色體結構信息和細胞間距離信息，進一步補充細胞的特征信息。實驗結果表明，scHiCSC對不同規模的數據集具有較強的泛化能力，性能優于現有的單細胞Hi-C聚類框架。雖然已經有針對單細胞Hi-C數據的線性或者非線性嵌入技術，但不同于基于單細胞RNA-seq數據的嵌入方法，目前尚未出現基于單細胞Hi-C數據的流形學習技術。單細胞Hi-C可以以矩陣、圖、向量的不同形式進行表達，這允許研究者們基于這些數據集進行更多的數據嵌入研究。
4 基于單細胞Hi-C數據的應用研究
對于器官/組織發育而言，構建人體細胞圖譜是了解器官功能的基礎[42]。單細胞組學在腦科學、腫瘤異質性、腫瘤免疫學、腫瘤耐藥和藥物開發等方面有一定的應用[43-44]。胚胎干細胞可以通過將自身逐漸暴露于特定因子，使用譜系特異性基因表達程序來產生功能細胞類型。Khateb等[45]使用單細胞Hi-C、RNA-seq、ATAC-seq等進行了胚胎干細胞多能性的染色質構像、單/體細胞基因表達、組蛋白修飾和可及性分析，獲得了前體細胞前中胚層的fate，并進一步進行了肌源性和神經源性分化，揭示了基因組和轉錄特征，并確定了指導初始Pax7表達和肌源性和神經源性過程激活的調控區域。Wang等[46]開發了EagleC工具，在慢性髓系白血病細胞系HAP1和K562上可預測染色體間移位的結構變異。也就是說，基于單細胞Hi-C數據，可進行基因結構變異的檢測，從而進行基因疾病的治療。
免疫細胞在分化和免疫應答過程中動態地重新編程其基因表達譜。然而，潛在的機制仍然難以捉摸。Lin等[47]開發了一種單細胞Hi-C方法，系統地描繪了巨噬細胞在這些過程中的三維基因組和動態表觀遺傳圖譜，并提出了“無序度”來測量拓撲關聯結構域內的基因組組織模式，這與單個細胞中的染色質表觀遺傳狀態、基因表達和染色質結構變異性相關。作者還發現NF-kappa B在結核分枝桿菌感染時啟動了系統性染色質構象重組。因此，單細胞Hi-C技術可以為研究細胞發育、細胞分化周期提供結構信息。
5 總結與展望
單細胞分辨率下分析三維基因組，為研究高階染色質結構和功能的異質性和動態性提供了獨特的機會。然而，單細胞Hi-C數據集的稀疏性、噪聲性和高維性等問題也帶來了新的分析挑戰。正如我們在這篇綜述中所討論的，計算方法已經專門用于單細胞Hi-C數據的分析，包括數據處理、降維、接觸圖插補、三維基因組結構建模和其他下游分析（例如，在不同尺度下調用三維基因組特征）。迄今為止，大多數單細胞Hi-C數據的計算方法都集中在降維和多尺度結構建模上，而其他類型的分析和解釋的重要方向需要更多的工作。
除了專門設計用于分析單細胞Hi-C數據的計算方法外，開發能夠將單細胞Hi-C數據與其他單細胞數據，例如單細胞RNA-seq數據、單細胞ATAC-seq數據等進行聯合分析的框架也至關重要。例如，聯合分析不同單細胞基因組數據的方法有望極大地促進對三維基因組結構、表觀基因組特征和轉錄活動的系統級理解。然而，單細胞Hi-C數據的極端稀疏性以及缺乏將單細胞Hi-C數據與其他單細胞分析相結合的系統評估，使得可靠的多模式分析具有相當的挑戰性。一個能同時分析單細胞三維基因組和其他表觀基因組以及功能特征的計算工具，將提供單細胞基因組結構和功能之間的直接關聯。總之，作為三維基因組和表觀基因組學研究中令人振奮的新方向，單細胞Hi-C和其他單細胞技術有望為核基因組結構和功能提供獨特和新穎的見解，這也將為開發新的計算方法創造機會。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：龔海燕查閱、整理相關文獻資料，麻付強起草論文，張曉彤負責論文整體框架和校對。
本文附件見本刊網站的電子版本（biomedeng.cn）。
