Hi-C 是一種基于測序的方法，用于分析全基因組染色質(zhì)互作。它已廣泛應(yīng)用于研究各種生物學(xué)問題，如基因調(diào)控、染色質(zhì)結(jié)構(gòu)、基因組組裝等。Hi-C 實驗涉及一系列生物化學(xué)反應(yīng)，可能會在輸出中引入噪聲。隨后的數(shù)據(jù)分析也會產(chǎn)生影響最終輸出噪聲：互作矩陣，其中矩陣中的每個元素表示基因組任意兩個區(qū)域之間的互作強度。因此，Hi-C 數(shù)據(jù)分析的關(guān)鍵步驟是消除此類噪聲，該步驟也稱為 Hi-C 數(shù)據(jù)歸一化。

歸一化方法

為了歸一化 Hi-C matrix，可以采用以下多種方法^[1]：

• Iterative Correction (IC) :該方法通過消除實驗過程中的偏差來歸一化原始接觸圖。這是一種矩陣平衡的方法，但是，在歸一化的情況下，行和列的總和不等于1。
• Knight-Ruiz Matrix Balancing (KR):Knight-Ruiz (KR) 矩陣平衡是一種歸一化對稱矩陣的快速算法。歸一化后獲得雙隨機矩陣。在這個矩陣中，行和列的總和等于一。
• Vanilla-Coverage (VC) :該方法首先用于染色體間圖譜。后來 Rao 等人，2014 年將其用于染色體內(nèi)圖譜。這是一種簡單的方法，首先將每個元素除以相應(yīng)行的總和，然后除以相應(yīng)列的總和。
• Median Contact Frequency Scaling (MCFS):此方法可用于使用兩個位置/坐標之間的特定距離的中值接觸值來歸一化接觸圖。首先，計算每個距離的中值距離接觸頻率。隨后，觀察到的接觸頻率除以根據(jù)兩個位置之間的距離獲得的中值接觸頻率。

方法詳解

早期的 Hi-C 數(shù)據(jù)歸一化方法主要關(guān)注引起噪聲的顯性因素。切割酶位點、基因組映射、GC 含量等因素使測序讀數(shù)在基因組中分布不均勻，從而在計算成對互作時引入偏差。根據(jù)這些想法，Imakaev 等人提出了一種能夠“implicitly”處理所有噪聲源的方法。背后的想法實際上非常簡單：因為 Hi-C 理論上是一個無偏的實驗，所有基因組區(qū)域的“visibility”應(yīng)該是相等的。另一個假設(shè)是 Hi-C 實驗的所有偏差都是一維且可分解的?；谶@些假設(shè)，一個解決方案是將原始互作矩陣分解為兩個一維偏差和一個行和列之和為相同值的歸一化矩陣的乘積。

Imakaev提出的方法在矩陣理論中也稱為矩陣平衡。矩陣平衡是一個古老的數(shù)學(xué)問題，可以追溯到一個世紀前。人們已經(jīng)開發(fā)了許多矩陣平衡方法，包括普通覆蓋率 (VC)、Sinkhorn & Knopp (SK) 和 Knight & Ruiz 方法 (KR)。

VC是通過將矩陣的每個元素除以其行和和列和來完成的，以去除每個位點的不同測序覆蓋度。 VC可以被認為是SK方法的單次迭代。在SK中，重復(fù)執(zhí)行VC過程，直到所有行和列的總和為相同的值。 S&K 這個過程會收斂。

Rao 等人回顧了所有矩陣平衡方法，并將 KR 方法引入 Hi-C 數(shù)據(jù)?；贙&R的原始論文，KR方法比SP快幾個數(shù)量級，這使得它適合平衡高分辨率矩陣。實際上，即使在 10kb 分辨率下，ICE 的 SP 實現(xiàn)也非?？臁８鶕?jù)我的經(jīng)驗，ICE 和 KR 之間的速度差異可以忽略不計。然而，KR 有一個缺點，即當矩陣太稀疏時，KR 過程可能無法收斂。在我的研究中，當我使用 Juicer tools 在低測序數(shù)據(jù)集上生成 KR 歸一化矩陣得到了一個空矩陣，這種情況發(fā)生了幾次。

矩陣平衡的算法其實并不難，我們?nèi)绾斡嬎?Hi-C 互作矩陣的平衡矩陣呢？下面的Python類中實現(xiàn)了VC和SP方法。對于小矩陣來說，這種實現(xiàn)速度很快。

class?HiCNorm(object):
????def?__init__(self,?matrix):
????????self.bias?=?None
????????self.matrix?=?matrix
????????self.norm_matrix?=?None
????????self._bias_var?=?None

????def?iterative_correction(self,?max_iter=50):
????????mat?=?np.array(self.matrix,?dtype=float)
????????row_sum?=?np.sum(mat,?axis=1)
????????low_count?=?np.quantile(row_sum,?0.15)
????????mask_row?=?row_sum?<?low_count
????????mat[mask_row,?:]?=?0
????????mat[:,?mask_row]?=?0

????????self.bias?=?np.ones(mat.shape[0])
????????self._bias_var?=?[]
????????x,?y?=?np.nonzero(mat)
????????mat_sum?=?np.sum(mat)

????????for?i?in?range(max_iter):
????????????bias?=?np.sum(mat,?axis=1)
????????????bias_mean?=?np.mean(bias[bias?>?0])
????????????bias_var?=?np.var(bias[bias?>?0])
????????????self._bias_var.append(bias_var)
????????????bias?=?bias?/?bias_mean
????????????bias[bias?==?0]?=?1

????????????mat[x,?y]?=?mat[x,?y]?/?(bias[x]*bias[y])
????????????new_sum?=?np.sum(mat)
????????????mat?=?mat?*?(mat_sum?/?new_sum)

????????????self.bias?=?self.bias?*?bias?*?np.sqrt(new_sum?/?mat_sum)

????????self.norm_matrix?=?np.array(self.matrix,?dtype=float)
????????self.norm_matrix[x,?y]?=?self.norm_matrix[x,?y]?/?(self.bias[x]?*?self.bias[y])

????def?vanilla_coverage(self):
????????self.iterative_correction(max_iter=1)

上述算法背后的想法是，我們首先將偏差設(shè)置為矩陣每行的總和，并將每個矩陣元素除以其行和列的偏差。重復(fù)這兩個步驟直到滿足收斂標準。我們可以使用偏差的方差（self.bias）來監(jiān)控平衡過程的收斂性（如下圖所示）。

原始互作矩陣、通過 SP 方法和 VC 方法歸一化的矩陣繪制為熱圖，如下所示。我們可以看到，歸一化矩陣中遠離對角線的區(qū)域比原始矩陣更干凈，但我們幾乎看不到 SP 和 VC 方法之間的差異。在實踐中，我們在歸一化之前預(yù)先過濾具有非常小的值的行。上面的腳本通過將這些行的元素設(shè)置為零來過濾掉總和低于所有行總和的 15 分位數(shù)的行。

然而，我們可以通過檢查相同距離的互作的相關(guān)性來量化 SP 和 VC 歸一化方法的差異。為此，我們提取并計算兩個矩陣的第 d 對角線的相關(guān)性，其中 d 是兩個基因組區(qū)域的距離（在 bin 處）。從下圖可以看出，雖然所有三種方法在長距離（>10 Mb）下都類似于原始矩陣，但 SP 與原始矩陣稍微相似。三種方法的成對比較表明，SP 和 VC 高度相似，只是迭代次數(shù)不同。 KR 比 VC 更類似于 SP，但差異可以忽略不計。 [注：KR矩陣是用straw API從.hic文件中獲得的]。

實際上，Hi-C 歸一化僅通過染色體內(nèi)、僅染色體間或兩者來完成。僅對于染色體內(nèi)，分別在每條染色體上進行 KR 或 ICE。僅對于染色體間，獲得全基因組矩陣并從中去除染色體內(nèi)互作。當包括染色體間相互作用時，高分辨率歸一化需要大量內(nèi)存。因此，歸一化通常在全基因組矩陣上以 25kb 或 50kb 分辨率進行。

一個自然的問題是，去除染色體間互作是否會對染色體內(nèi)互作的產(chǎn)生影響？為了回答這個問題，我對所有互作點和僅染色體內(nèi)互作點進行了 SP 標準化。同樣，通過所有互作歸一化和僅通過染色體內(nèi)互作歸一化之間的差異非常小。 Juicer tools 和 Cooler 都默認使用所有觸點進行歸一化。

總結(jié)

那么我應(yīng)該對數(shù)據(jù)使用哪種歸一化方法？

答案可能是這真的不重要。正如 Rao 等人和我們的分析所示，在 VC、ICE(SP)、KR 和其他幾種矩陣平衡方法之間觀察到高度相關(guān)性。 Rao等人在他們的研究中進一步指出，循環(huán)調(diào)用不受不同歸一化方法的影響。他們甚至表明，當對原始數(shù)據(jù)調(diào)用峰值時，循環(huán)幾乎相同，這讓我懷疑我們是否需要矩陣歸一化。雖然Rao等人確實提到VC傾向于過度修正互作圖并提出VCSQ來補充VC，但我還沒有看到有人拿出明確的證據(jù)表明KR/ICE優(yōu)于VC/VCSQ。每種方法的優(yōu)缺點總結(jié)如下：

總之，仍然建議將歸一化作為實踐的一部分。選擇 KR 還是 ICE 在很大程度上取決于您使用的流程（Juicer tools 或 Cooler），但如果您的測序深度相對較低，推薦 ICE，因為它可以保證收斂。