基因組中的大部分 DNA 以染色體的形式緊密盤繞在細(xì)胞核內(nèi)，基因進(jìn)行復(fù)制和轉(zhuǎn)錄時，高度折疊的染色質(zhì)結(jié)構(gòu)需要暴露出 DNA 序列，這段暴露的區(qū)域就是染色質(zhì)開放區(qū)域（也叫開放染色質(zhì) open chromatin），這個區(qū)域可以供轉(zhuǎn)錄因子和其他調(diào)控元件結(jié)合。暴露的染色質(zhì)開放區(qū)域會包含啟動子、增強子、絕緣子、沉默子等順式調(diào)控元件，其可以被反式作用因子接近結(jié)合的特性，為染色質(zhì)的可及性（也叫染色質(zhì)開放性 chromatin accessibility ）。
DNA 開放區(qū)域與基因表達(dá)類似，在不同的組織細(xì)胞中、在發(fā)育或細(xì)胞的不同時期動態(tài)變化，呈現(xiàn)時空變化特異性，通過 ATAC-seq（Assay for Transposase-Accessible Chromatin using sequencing）技術(shù)可實現(xiàn) bulk 樣本染色質(zhì)可及性的分析從而進(jìn)行相關(guān)的表觀研究。
隨著近幾年單細(xì)胞水平技術(shù)平臺的發(fā)展，采用 10x Genomics Chromium Single Cell ATAC 等產(chǎn)品平臺可以得到單個細(xì)胞中染色質(zhì)開放區(qū)域圖譜數(shù)據(jù)，也可采用 10x Genomics Single Cell Multiome ATAC + Gene Expression 類似的多組學(xué)平臺同時獲得單細(xì)胞水平的基因表達(dá)數(shù)據(jù)和表觀調(diào)控數(shù)據(jù)，從而進(jìn)行特定細(xì)胞類型下更加精細(xì)化基因表達(dá)的調(diào)控機制研究。

今天為大家分享 scATAC 單個平臺數(shù)據(jù)的分析內(nèi)容和分析思路，希望對各位老師有所幫助！

scATAC 分析思路

1、質(zhì)控

TileMatrix：fragments 可以理解為 scATAC 測序得到 DNA 開放區(qū)域片段，根據(jù) fragments 的基因組定位信息，如以 5kb 為一個滑動窗口 bin，依次統(tǒng)計單個細(xì)胞在每個 bin 中 fragments 的數(shù)目，最終形成 TileMatrix，詳見 ArchR分 析軟件[1]。

PeakMatrix：對全部細(xì)胞的 fragments 進(jìn)行 peak calling，過濾掉低可信度的 DNA 開放區(qū)域數(shù)據(jù)，保留在多個細(xì)胞中保守的染色質(zhì)開放區(qū)域，即為 peak，最終生成 Peak-Barcode 可及性矩陣。大部分文獻(xiàn)和軟件（如 cellranger-atac[2]?和 Signac[3]）都是默認(rèn)基于 PeakMatirx 進(jìn)行分析，后續(xù)觀察 marker gene 上下游區(qū)域尤其是 TSS 區(qū)域在不同細(xì)胞群間的分布差異，指導(dǎo)細(xì)胞類型定義。單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)的基因 feature 數(shù)目一般在2~3萬，而 scATAC 數(shù)據(jù) peak calling 結(jié)果一般會包含~10萬左右的 peaks。

圖中淡黃色部分有有效的細(xì)胞 barcodes，藍(lán)色的點表示非細(xì)胞 barcodes，橫坐標(biāo)表示每個 barcode fragments 的數(shù)目，縱坐標(biāo)為單個 barcode 內(nèi)屬于 peak 部分的 fragments 占全部 fragments 的比例。從上圖可以看出 cellranger-atac 算法認(rèn)為只有包含 fragmengts 要超過一定的閾值，且保守的可信 peaks 區(qū)的 fragments 的比例大于一定的閾值的 barcodes 才是有效的細(xì)胞。

TSSenrichment score：?也是文獻(xiàn)中常見的一個細(xì)胞識別過濾參數(shù)，細(xì)胞保持活性需要一定數(shù)量基因的表達(dá)翻譯，而基因轉(zhuǎn)錄時就需要開放 TSS 區(qū)域，好結(jié)合轉(zhuǎn)錄因子、DNA 聚合酶等轉(zhuǎn)錄起始復(fù)合物，一般情況下 TSS 區(qū)域的 fragments 能占到細(xì)胞全部 fragments 的25%~35%，如上圖所示，一般要求有效細(xì)胞 barcode 的 TSS 富集得分>4且 barcode 包含 fragments>1000，也可設(shè)置更為嚴(yán)格的過濾參數(shù)，如 TSS 富集得分 >6 或者 >8。

scATAC 數(shù)據(jù)過濾的其它參數(shù)有：

nucleosome_signal：反映 DNA 開放區(qū)域的長度信息，DNA 開放區(qū)域的長度通常會在 1~2 個核小體長度（147bp）范圍內(nèi)，所以有些文獻(xiàn)會過濾掉 nucleosome_signal>4的 barcodes。

blacklist regions ratio：ENCODE 數(shù)據(jù)庫中收集了 blacklist 區(qū)域（基因組總反常的或者無論在二代測序的哪個實驗中都是高信號的區(qū)域，排除掉這些區(qū)域可以更好的分析功能基因組數(shù)據(jù)），因此要求 barcode 中比對到 blacklist 的 fragments 的比例要低于一定的閾值。

與單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)雙細(xì)胞過濾分析類似，如果一個 barcode 包含的 fragments 太多，則此 barcode 可能包含了兩個或多個細(xì)胞，除了設(shè)定 fragments 數(shù)目閾值進(jìn)行簡單過濾外，同樣可采用 demuxlet[7]?進(jìn)行雙細(xì)胞過濾，注意（圖3）中 Doublet Enrichment 得分并不等同于細(xì)胞數(shù)目，僅按照此得分從大到小排序，然后按比例進(jìn)行過濾。

2、轉(zhuǎn)錄因子可及性分析

TF-motif Matrix：一個轉(zhuǎn)錄因子（Transcription Factor）通常包含一個可結(jié)合 DNA 基序（motif），少數(shù)情況下一個轉(zhuǎn)錄因子會有兩個或多個 motif，通過對近 10 萬 peaks 的反向搜索，統(tǒng)計每個 TF-motif 在每個細(xì)胞中可以結(jié)合的 fragments 數(shù)目，從而得到單細(xì)胞層面轉(zhuǎn)錄因子水平的可及性矩陣 TF-motifMatrix，常見的分析軟件有 chromVAR[4]，常用的JASPAR[5]?核心庫包含的 TF~motif 在 600 左右，因此從 PeakMatrix 轉(zhuǎn)化為 TF-motif Matrix 也可以看作為一個數(shù)據(jù)降維的過程。

注：縱坐標(biāo)表示每個 TF motif 在全部細(xì)胞中的可及性變化程度，橫坐標(biāo)為排秩后的轉(zhuǎn)錄因子。

3、細(xì)胞聚類

除了常規(guī)的基于 TileMatrix、PeakMatrix 進(jìn)行降維，細(xì)胞聚類分群外， Jason D. Buenrostro 2018 年在?Cell?發(fā)表的關(guān)于人類造血細(xì)胞分化中表觀調(diào)控的文獻(xiàn)，采用了 TF motif Matirx 數(shù)據(jù)進(jìn)行細(xì)胞聚類分析（如下圖所示），當(dāng)然基于 GeneScore Matirx 參考單細(xì)胞轉(zhuǎn)錄分析流程進(jìn)行聚類分析也是一種思路。

初始聚類->peak calling->再聚類模式：上圖為 Ansuman T. Satpathy 在?NATURE BIOTECHNOLOGY?發(fā)表的關(guān)于人類免疫細(xì)胞發(fā)育的表觀圖譜文獻(xiàn)的分析思路。簡介如下：1. 保留 Unique fragments >1000且 TSS 區(qū)域富集得分>8的 barcode 為有效細(xì)胞；2. 基因組位置上每 2.5kb 為一個 bin，生成 TileMatrix，選擇變化最大最高的 top2 萬 bins 進(jìn)行初始聚類；3. 依據(jù)初始分類結(jié)果，每 200 個細(xì)胞合并成一個 pseudo bulk ATAC 數(shù)據(jù)；4. 對每個 Cluster 的 pseudo bulk ATAC 數(shù)據(jù)分別 call peaks，最后合并全部的 peaks 數(shù)據(jù)；5. 生成最終的 PeakMatrix、降維、聚類、細(xì)胞類型注釋及后續(xù)分析。
采用類似思路三輪聚類作為最終結(jié)果的還有 Alexandro E. Trevino 2020 年在?Cell?上發(fā)表的關(guān)于人類大腦皮層發(fā)育的表觀基因表達(dá)調(diào)控研究的文獻(xiàn)。此分析思路總結(jié)如下：一方面，合并單細(xì)胞數(shù)據(jù)生成 pseudo bulk ATAC 數(shù)據(jù)，可以有效減少由于單個細(xì)胞 ATAC 數(shù)據(jù)捕獲不足造成的數(shù)據(jù)缺失問題；另一方面，每個初始 cluster 重新 call peak，可以得到細(xì)胞群間特異的 peaks，便于發(fā)現(xiàn)細(xì)胞特異的表觀調(diào)控機制。

4、細(xì)胞類型注釋

Marker 基因方法：在 scATAC 細(xì)胞分群結(jié)果的基礎(chǔ)上，進(jìn)行類間差異 peak 可及性分析，根據(jù) marker 基因 TSS 區(qū)域或整個基因區(qū)域的可及性類間差異，對細(xì)胞群進(jìn)行人工注釋。

基因活性方法：GeneScoreMatrix?基因活性得分矩陣，考慮 peak 與 TSS/genebody 的距離、測序深度等因素采用廣義線性模型將 PeakMatrix 矩陣轉(zhuǎn)換為 GeneScore Matrix，常用的軟件有 Cicero[6]、Signac[3]?和 ArchR[1]?；蚧钚缘梅挚梢岳斫鉃榛诒碛^數(shù)據(jù)預(yù)測出的基因表達(dá)值，可以通過觀察 markergene 在細(xì)胞群中活性得分值的特異分布，進(jìn)行細(xì)胞類型的人工注釋。也可采用 SingleR 和 Cellassign 等軟件進(jìn)行細(xì)胞類型注釋。

聯(lián)合 scRNA 轉(zhuǎn)移類標(biāo)簽方法：根據(jù)實驗設(shè)計尋找最接近且已經(jīng)注釋好細(xì)胞類型的 scRNA 數(shù)據(jù)，基于 Gene Score 數(shù)據(jù)，采用 Seurat TransferData 的方法，獲得 scATAC 的細(xì)胞類型注釋。

5、共可及性分析

peak-peak 共可及性分析，常用來研究基因組上的順勢作用元件間（如啟動子和增強子）的相互作用，可采用 Cicero[6]包實現(xiàn)。

如上圖所示，橫坐標(biāo)是預(yù)測到的共可及性網(wǎng)絡(luò)區(qū)間內(nèi)的 peaks 和注釋到的基因，縱坐標(biāo)是共可及性值，圖中 0.2 處的虛線是存在共可及性關(guān)系的閾值線，紫色弧線是 peak 之間的相互作用連接線。

6、擬時間分析

基于基因活性得分矩陣，采用 Monocle[8]?軟件可實現(xiàn) scATAC 數(shù)據(jù)的擬時間分析（又稱趨勢分析）。下圖為 scATAC 的細(xì)胞在 Monocle 降維結(jié)果上 pseudotime 的分布圖。

參考文獻(xiàn)

1. Granja JM, Corces MR et al., ArchR is a scalable software package for integrative single-cell chromatin accessibility analysis.?Nature Genetics?(2021).

2. https://support.10xgenomics.com/single-cell-atac/software/pipelines/latest/algorithms/overview.

3. Tim Stuart, Avi Srivastava, et al. Single-cell chromatin state analysis with Signac.?Nat Methods?. 2021 Nov;18(11):1333-1341. HelenaTodorov. et al.

4. Schep AN, Wu B, Buenrostro JD and Greenleaf WJ. chromVAR: inferring transcription-factor-associated accessibility from single-cell epigenomic data.?Nature Methods. doi: 10.1038/nmeth.4401.

5. Khan, A. JASPAR 2018: update of the open-access database of transcription factor binding profiles and its web framework.?Nucleic Acids Res. 2018; 46: D260–D266.

6. Hannah A. Pliner, Jonathan S. Packer, et al. Cicero Predicts cis-Regulatory DNA Interactions from Single-Cell Chromatin Accessibility Data.?Molecular Cell. TECHNOLOGY| VOLUME 71, ISSUE 5, P858-871. E8, SEPTEMBER 06, 2018

7. Kang et al.,Multiplexed droplet single-cell RNA-sequencing using natural genetic variation.?Nature Biotechnology?2017.

8. Cole Trapnell, Davide Cacchiarelli. The dynamics and regulators of cell fate decisions are revealed by pseudotemporal ordering of single cells.?Nature Biotechnology?2014.