scATAC分析神器ArchR初探-簡(jiǎn)介（1）
scATAC分析神器ArchR初探-ArchR進(jìn)行doublet處理（2）
scATAC分析神器ArchR初探-創(chuàng)建ArchRProject（3）
scATAC分析神器ArchR初探-使用ArchR降維（4）
scATAC分析神器ArchR初探--使用ArchR進(jìn)行聚類（5）
scATAC分析神器ArchR初探-單細(xì)胞嵌入（6）
scATAC分析神器ArchR初探-使用ArchR計(jì)算基因活性值和標(biāo)記基因（7）
scATAC分析神器ArchR初探-scRNA-seq確定細(xì)胞類型（8）
scATAC分析神器ArchR初探-ArchR中的偽批次重復(fù)處理（9）
scATAC分析神器ArchR初探-使用ArchR-peak-calling（10）
scATAC分析神器ArchR初探-使用ArchR識(shí)別標(biāo)記峰（11）
scATAC分析神器ArchR初探-使用ArchR進(jìn)行主題和功能豐富(12)
scATAC分析神器ArchR初探-利用ArchR豐富ChromVAR偏差(13)
scATAC分析神器ArchR初探-使用ArchR進(jìn)行足跡（14）
scATAC分析神器ArchR初探-使用ArchR進(jìn)行整合分析（15）
scATAC分析神器ArchR初探-使用ArchR進(jìn)行軌跡分析（16）

簡(jiǎn)介

近來單細(xì)胞技術(shù)發(fā)展迅速，不僅scRNA技術(shù)應(yīng)用逐漸成熟，scATAC技術(shù)也在開始規(guī)模應(yīng)用，分析工具也多種多樣阀溶，這里介紹由ATAC發(fā)明者greenleaf實(shí)驗(yàn)室開發(fā)的單細(xì)胞ATAC分析工具ArchR使用說明，ArchR：對(duì)單細(xì)胞染色質(zhì)可及性數(shù)據(jù)的可靠且可擴(kuò)展的分析杈抢，ArchR是功能齊全的軟件套件京腥，用于分析單細(xì)胞染色質(zhì)可訪問性數(shù)據(jù)蜜另。它設(shè)計(jì)用于處理成千上萬個(gè)單細(xì)胞碉纳，而沒有大的內(nèi)存或計(jì)算需求勿负，與10x Genomics Chromium系統(tǒng)等商業(yè)平臺(tái)可達(dá)到的實(shí)驗(yàn)規(guī)模保持同步。這部分介紹如何將數(shù)據(jù)導(dǎo)入ArchR劳曹，以及如何創(chuàng)建ArrowFiles（ArchR分析的基本單元）

1.1關(guān)于ATAC-seq術(shù)語的簡(jiǎn)要入門

任何ATAC-seq實(shí)驗(yàn)的最基本組成部分是“ 片段”奴愉。在ATAC-seq中，片段是指通過兩次轉(zhuǎn)座事件產(chǎn)生的可測(cè)序DNA分子铁孵。使用配對(duì)末端測(cè)序?qū)ζ蔚拿總€(gè)末端進(jìn)行測(cè)序锭硼。基于Tn5的插入偏移量調(diào)整推斷的片段起點(diǎn)和終點(diǎn)的單堿基位置蜕劝。如先前報(bào)道檀头，Tn5轉(zhuǎn)座酶以同源二聚體的形式與DNA結(jié)合，兩個(gè)Tn5分子之間具有9 bp的DNA岖沛。因此暑始，每個(gè)Tn5同型二聚體結(jié)合事件都會(huì)產(chǎn)生兩個(gè)插入，間隔為9 bp烫止。因此，“可及”位點(diǎn)的實(shí)際中心點(diǎn)位于Tn5二聚體的正中心戳稽，而不是每個(gè)Tn5插入的位置馆蠕。為了解決這個(gè)問題，我們將偏移量應(yīng)用于單個(gè)Tn5插入惊奇，將正鏈插入事件調(diào)整為+4 bp互躬，將負(fù)鏈插入事件調(diào)整為-5 bp。ATAC-seq的原始描述颂郎。因此吼渡，在ArchR中，“ 片段 ”是指一個(gè)表或基因組范圍對(duì)象乓序，其中包含染色體寺酪，偏移量調(diào)整后的單堿基染色體起始位置，偏移量調(diào)整后的單堿基染色體終點(diǎn)位置以及與每個(gè)序列化片段相對(duì)應(yīng)的唯一細(xì)胞條形碼ID 替劈。類似地寄雀，“ 插入 ”是指在可到達(dá)部位的正中央的偏移調(diào)整后的單堿基位置。

1.2為什么要使用ArchR陨献？

那里有多種用于單細(xì)胞ATAC序列分析的工具盒犹，那么為什么要使用ArchR？最重要的是，ArchR提供功能并支持其他工具無法做到的分析：

此外急膀，由于對(duì)數(shù)據(jù)結(jié)構(gòu)和并行化方法進(jìn)行了重大優(yōu)化沮协，構(gòu)成了ArchR軟件的基礎(chǔ)，因此與其他可用工具相比卓嫂，ArchR更快并且使用的內(nèi)存更少慷暂。當(dāng)分析超過70,000個(gè)單元時(shí)，某些工具需要高性能的計(jì)算環(huán)境命黔，超過128 GB的可用內(nèi)存（OoM =內(nèi)存不足）呜呐。

1.3什么是箭頭文件/ ArchRProject？

ArchR中分析項(xiàng)目的基本單位稱為箭頭文件背桐。每個(gè)Arrow文件都存儲(chǔ)與單個(gè)樣本關(guān)聯(lián)的所有數(shù)據(jù)（即优烧，元數(shù)據(jù)，可訪問的片段和數(shù)據(jù)矩陣）链峭。在此畦娄，“單個(gè)樣本”將是所需的最詳細(xì)的分析單位（例如，特定條件的單個(gè)重復(fù)）弊仪。在創(chuàng)建過程中以及執(zhí)行附加分析時(shí)熙卡，ArchR更新和編輯每個(gè)Arrow文件以包含附加信息層。值得注意的是励饵，對(duì)于ArchR再膳，Arrow文件實(shí)際上只是磁盤上存儲(chǔ)的外部文件的路徑。更明確地說曲横，箭頭文件不是存儲(chǔ)在內(nèi)存中的R語言對(duì)象喂柒，而是存儲(chǔ)在磁盤上的HDF5格式的文件不瓶。因此，我們使用ArchRProject對(duì)象灾杰，將這些Arrow文件關(guān)聯(lián)在一起蚊丐，成為一個(gè)可以在R中快速訪問的單個(gè)分析框架。此ArchRProject對(duì)象尺寸很小艳吠，并存儲(chǔ)在內(nèi)存中麦备。

某些操作可以直接在Arrow文件上執(zhí)行，而其他操作則可以在上執(zhí)行ArchRProject昭娩，從而依次更新每個(gè)關(guān)聯(lián)的Arrow文件凛篙。由于Arrow文件存儲(chǔ)為大型HDF5格式文件，因此ArchR中的“ get-er”功能通過與ArchRProjectwhile“ add-er”功能進(jìn)行交互來檢索數(shù)據(jù)栏渺，或者（i）將數(shù)據(jù)直接添加到Arrow文件呛梆，（ii）直接添加數(shù)據(jù)到ArchRProject，或（iii）通過與交互數(shù)據(jù)添加到箭頭文件ArchRProject磕诊。

1.4 ArchR中的輸入文件類型

ArchR可以利用scATAC-seq數(shù)據(jù)的多種輸入格式填物，這種格式最常見的是片段文件和BAM文件的格式。片段文件是由tabix排序的文本文件霎终，包含每個(gè)scATAC-seq片段和相應(yīng)的單元ID滞磺，每行一個(gè)片段。BAM文件是按Tabix排序的二進(jìn)制文件莱褒，包含每個(gè)scATAC-seq片段击困，原始序列，蜂窩條形碼ID和其他信息广凸。使用的輸入格式將取決于使用的預(yù)處理管道阅茶。例如，10x Genomics Cell Ranger軟件返回片段文件炮障，而sci-ATAC-seq應(yīng)用程序?qū)⑹褂肂AM文件目派。ArchR使用“ scanTabix”讀取片段文件坤候，使用“ scanBam”讀取BAM文件胁赢。在此輸入過程中，對(duì)輸入進(jìn)行分塊白筹，然后將每個(gè)輸入塊轉(zhuǎn)換為基于壓縮表的片段表示形式智末，其中包含每個(gè)片段染色體，偏移校正后的染色體起始位置徒河，偏移校正后的染色體終止位置和細(xì)胞條形碼ID系馆。然后，將這些逐塊片段存儲(chǔ)在HDF5格式的臨時(shí)文件中顽照，以保留內(nèi)存使用率由蘑，同時(shí)保持對(duì)每個(gè)塊的快速訪問闽寡。最后，與每個(gè)染色體相關(guān)的所有塊均被讀取尼酿，組織并重新寫入一個(gè)稱為“片段”的HDF5組內(nèi)的Arrow文件爷狈。此預(yù)分塊過程使ArchR可以有效地處理非常大的輸入文件，并且內(nèi)存使用率較低裳擎，從而可以充分利用并行處理涎永。與每個(gè)染色體相關(guān)的所有塊均被讀取，組織并重寫為一個(gè)稱為“片段”的HDF5組中的Arrow文件鹿响。此預(yù)分塊過程使ArchR可以有效地處理非常大的輸入文件羡微，并且內(nèi)存使用率較低，從而可以充分利用并行處理惶我。與每個(gè)染色體相關(guān)的所有塊均被讀取妈倔，組織并重寫為一個(gè)稱為“片段”的HDF5組中的Arrow文件。此預(yù)分塊過程使ArchR可以有效地處理非常大的輸入文件指孤，并且內(nèi)存使用率較低启涯，從而可以充分利用并行處理。

1.5設(shè)置

我們要做的第一件事是更改到所需的工作目錄恃轩，設(shè)置我們要使用的線程數(shù)结洼，并加載我們的基因和基因組注釋。根據(jù)您本地環(huán)境的配置叉跛，您可能需要在中修改threads以下使用的數(shù)量addArchRThreads()松忍。默認(rèn)情況下，ArchR使用threads可用總數(shù)的一半筷厘，但是您可以根據(jù)需要手動(dòng)進(jìn)行調(diào)整鸣峭。如果使用Windows，則threads因?yàn)锳rchR中的并行處理是針對(duì)基于Unix的操作系統(tǒng)而構(gòu)建的酥艳，所以可用將自動(dòng)設(shè)置為1摊溶。

首先，我們加載ArchR庫充石。如果失敗莫换，則說明您尚未正確安裝ArchR，應(yīng)重新訪問安裝說明骤铃。

library(ArchR)

接下來拉岁，我們?cè)O(shè)置ArchR函數(shù)的默認(rèn)線程數(shù)。這是您在每個(gè)新的R會(huì)話期間都必須執(zhí)行的操作惰爬。我們建議將threads可用內(nèi)核總數(shù)設(shè)置為1/2到3/4喊暖。ArchR中的內(nèi)存使用量通常會(huì)隨所用線程的數(shù)量而定，因此允許ArchR使用更多線程也將導(dǎo)致更高的內(nèi)存使用量撕瞧。

addArchRThreads(threads = 16) 

將并行線程的默認(rèn)數(shù)量設(shè)置為16陵叽。

然后狞尔，我們?cè)O(shè)置用于基因和基因組注釋的基因組。如上所述巩掺，這是您在每個(gè)新的R會(huì)話期間都必須執(zhí)行的操作沪么。當(dāng)然，該基因組版本必須與用于比對(duì)的基因組版本匹配锌半。對(duì)于本教程中使用的數(shù)據(jù)禽车，我們將使用hg19參考基因組，但是ArchR本機(jī)支持其他基因組注釋和自定義基因組注釋刊殉，如下一節(jié)所述殉摔。

addArchRGenome("hg19")

將默認(rèn)基因組設(shè)置為Hg19。##設(shè)置基因組和基因注釋

ArchR需要基因和基因組注釋來執(zhí)行諸如計(jì)算TSS富集得分记焊，核苷酸含量和基因活性得分之類的事情逸月。由于我們的教程數(shù)據(jù)集使用已經(jīng)與hg19參考基因組對(duì)齊的scATAC-seq數(shù)據(jù)，因此在上一節(jié)中我們將“ hg19”設(shè)置為默認(rèn)基因組遍膜。但是碗硬，ArchR本地支持“ hg19”，“ hg38”瓢颅，“ mm9”和“ mm10”恩尾，您可以使用createGeneAnnotation()和createGenomeAnnotation()函數(shù)創(chuàng)建自己的基因組和基因注釋。

通過該addArchRGenome()功能可以簡(jiǎn)化向ArchR提供此信息的過程挽懦。此函數(shù)告訴ArchR翰意，對(duì)于當(dāng)前會(huì)話中的所有分析，它都應(yīng)使用genomeAnnotation和geneAnnotation與define 關(guān)聯(lián)ArchRGenome信柿。每個(gè)由本地支持的基因組均由BSgenome定義每個(gè)染色體的基因組坐標(biāo)和序列的GRanges對(duì)象冀偶，包含一組列入黑名單的區(qū)域的TxDb對(duì)象，定義所有基因的位置和結(jié)構(gòu)的OrgDb對(duì)象以及提供中心基因的對(duì)象組成基因標(biāo)識(shí)符渔嚷，并包含此標(biāo)識(shí)符與其他種類的標(biāo)識(shí)符之間的映射进鸠。

以下是如何為本地支持的基因組加載基因和基因組注釋的示例，以及有關(guān)其BSgenome和黑名單組件的信息形病。

在預(yù)編譯的版本hg19基因組ArchR使用BSgenome.Hsapiens.UCSC.hg19客年，TxDb.Hsapiens.UCSC.hg19.knownGene，org.Hs.eg.db窒朋，和用合并黑名單ArchR::mergeGR()從hg19 V2黑名單區(qū)和顯示高mappability到hg19核基因組線粒體地區(qū)從迦勒Lareau和Jason Buenrostro搀罢。要將全局基因組默認(rèn)設(shè)置為預(yù)編譯的hg19基因組蝗岖，請(qǐng)執(zhí)行以下操作：

addArchRGenome("hg19")

將默認(rèn)基因組設(shè)置為Hg19侥猩。

在預(yù)編譯的版本hg38基因組ArchR使用BSgenome.Hsapiens.UCSC.hg38，TxDb.Hsapiens.UCSC.hg38.knownGene抵赢，org.Hs.eg.db欺劳，和用合并黑名單ArchR::mergeGR()從hg38 V2黑名單區(qū)和顯示高mappability到hg38核基因組線粒體地區(qū)從迦勒Lareau和Jason Buenrostro唧取。要將全局基因組默認(rèn)設(shè)置為預(yù)編譯的hg38基因組，請(qǐng)執(zhí)行以下操作：

addArchRGenome("hg38")

將默認(rèn)基因組設(shè)置為Hg38划提。

在預(yù)編譯的版本MM9在ArchR用途的基因組BSgenome.Mmusculus.UCSC.mm9枫弟，TxDb.Mmusculus.UCSC.mm9.knownGene，org.Mm.eg.db鹏往，并使用合并淡诗，這是一個(gè)黑名單ArchR::mergeGR()從MM9 V1黑名單地區(qū)從Anshul Kundaje和顯示高mappability至MM9核基因組線粒體地區(qū)從迦勒Lareau和Jason Buenrostro。要將全局基因組默認(rèn)設(shè)置為預(yù)編譯的mm9基因組伊履，請(qǐng)執(zhí)行以下操作：

addArchRGenome("mm9")

將默認(rèn)基因組設(shè)置為Mm9韩容。

在預(yù)編譯的版本MM10基因組ArchR使用BSgenome.Mmusculus.UCSC.mm10，TxDb.Mmusculus.UCSC.mm10.knownGene唐瀑，org.Mm.eg.db群凶，和用合并黑名單ArchR::mergeGR()從MM10 V2黑名單區(qū)和顯示高mappability到MM10核基因組線粒體地區(qū)從迦勒Lareau和Jason Buenrostro。要將全局基因組默認(rèn)設(shè)置為預(yù)編譯的mm10基因組哄辣，請(qǐng)執(zhí)行以下操作：

addArchRGenome("mm10")

將默認(rèn)基因組設(shè)置為Mm10请梢。

1.5.1創(chuàng)建自定義ArchRGenome

如上所述，ArchRGenome由基因組注釋和基因注釋組成力穗。

要?jiǎng)?chuàng)建自定義基因組注釋毅弧，我們使用createGenomeAnnotation()。為此当窗，您將需要以下信息：

1. 甲BSgenome其中包含一個(gè)基因組序列信息的對(duì)象形真。這些通常是Bioconductor軟件包（例如，BSgenome.Hsapiens.UCSC.hg38）超全，可通過google輕松找到咆霜。
2. 一個(gè)GRanges基因組范圍對(duì)象，其中包含一組列入黑名單的區(qū)域嘶朱，這些區(qū)域?qū)⒂糜趶南掠畏治鲋羞^濾掉不需要的區(qū)域蛾坯。這不是必需的，但建議這樣做疏遏。有關(guān)如何創(chuàng)建黑名單的信息脉课，請(qǐng)參閱ENCODE黑名單上的出版物。

例如财异，如果我們想從果蠅中創(chuàng)建自定義基因組注釋倘零，我們將首先識(shí)別并安裝并加載相關(guān)BSgenome對(duì)象。

if (!requireNamespace("BSgenome.Dmelanogaster.UCSC.dm6", quietly = TRUE)){
  BiocManager::install("BSgenome.Dmelanogaster.UCSC.dm6")
}
library(BSgenome.Dmelanogaster.UCSC.dm6)

加載所需的軟件包：BSgenome
加載所需的軟件包：Biostrings
加載所需的軟件包：XVector
附加軟件包：'Biostrings'
以下對(duì)象從'package：base'中屏蔽：strsplit

然后戳寸，我們從該BSgenome對(duì)象創(chuàng)建基因組注釋呈驶。

genomeAnnotation <- createGenomeAnnotation(genome = BSgenome.Dmelanogaster.UCSC.dm6)

正在獲取基因組..
正在獲取chromSizes ..
正在獲取黑名單..
未下載黑名單！如果沒有繼續(xù)疫鹊，請(qǐng)注意下游的偏差袖瞻。

檢查該對(duì)象顯示了ArchR基因組注釋的組成部分司致。

genomeAnnotation

長(zhǎng)度為3的列表
名稱（3）：基因組chromSizes黑名單

要?jiǎng)?chuàng)建自定義基因注釋，我們使用createGeneAnnotation()聋迎。為此脂矫，您將需要以下信息：

1. 甲TxDb從Bioconductor的對(duì)象（轉(zhuǎn)錄數(shù)據(jù)庫），其包含用于基因/轉(zhuǎn)錄物的坐標(biāo)的信息霉晕。
2. OrgDb來自Bioconductor 的對(duì)象（生物數(shù)據(jù)庫）庭再，提供了一個(gè)統(tǒng)一的框架來映射基因名稱和各種基因標(biāo)識(shí)符。

繼續(xù)以果蠅（Drosophila melanogaster）為例牺堰，我們首先安裝并加載相關(guān)對(duì)象TxDb和OrgDb對(duì)象佩微。

if (!requireNamespace("TxDb.Dmelanogaster.UCSC.dm6.ensGene", quietly = TRUE)){
  BiocManager::install("TxDb.Dmelanogaster.UCSC.dm6.ensGene")
}
if (!requireNamespace("org.Dm.eg.db", quietly = TRUE)){
  BiocManager::install("org.Dm.eg.db")
}
library(TxDb.Dmelanogaster.UCSC.dm6.ensGene)

加載所需的軟件包：GenomicFeatures
加載所需的軟件包：AnnotationDbi

library(org.Dm.eg.db)

然后我們創(chuàng)建基因注釋對(duì)象。

geneAnnotation <- createGeneAnnotation(TxDb = TxDb.Dmelanogaster.UCSC.dm6.ensGene, OrgDb = org.Dm.eg.db)

正在獲取基因..
確定的注釋樣式= ENSEMBL
正在獲取外顯子..
正在獲取TSS ..

檢查該對(duì)象顯示了ArchR基因注釋的組成部分萌焰。

geneAnnotation

長(zhǎng)度為3的列表
名稱（3）：基因外顯子TSS

另外哺眯，如果您沒有TxDb和OrgDb對(duì)象，則可以geneAnnotation根據(jù)以下信息創(chuàng)建一個(gè)對(duì)象：

1. “基因”對(duì)象- GRanges包含基因坐標(biāo)的對(duì)象（開始到結(jié)束）扒俯。該對(duì)象的符號(hào)列必須與下面描述的“外顯子”對(duì)象的符號(hào)列匹配奶卓。
2. “外顯子”對(duì)象- GRanges包含基因外顯子坐標(biāo)的對(duì)象。必須具有與上述“基因”對(duì)象的符號(hào)列匹配的符號(hào)列撼玄。
3. 一個(gè)GRanges包含standed轉(zhuǎn)錄起始位點(diǎn)（對(duì)象TSS）的坐標(biāo)夺姑。

geneAnnotation <- createGeneAnnotation(
  TSS = geneAnnotation$TSS, 
  exons = geneAnnotation$exons, 
  genes = geneAnnotation$genes
)

1.5.2在ArchR中使用非標(biāo)準(zhǔn)基因組

ArchR實(shí)施了一些檢查，以防止您執(zhí)行我們認(rèn)為“超出正常范圍”的事情掌猛。這些檢查之一會(huì)強(qiáng)制seqnames您的基因組注釋以“ chr”開頭盏浙。在大多數(shù)應(yīng)用中都是如此，但是有些基因組（例如玉米）不使用“ chr”作為不同染色體的前綴荔茬。要對(duì)此類基因組進(jìn)行任何ArchR分析废膘，您必須告訴ArchR忽略染色體前綴。為此慕蔚，您必須addArchRChrPrefix(chrPrefix = FALSE)在創(chuàng)建Arrow文件之前運(yùn)行丐黄。這將在當(dāng)前R會(huì)話中全局關(guān)閉對(duì)seqnames的“ chr”前綴的要求。請(qǐng)記住孔飒，ArchR默認(rèn)將染色體/序列名稱轉(zhuǎn)換為字符灌闺。因此，如果你seqnames僅僅是整數(shù)坏瞄，你需要提供他們的角色桂对，即 useSeqnames = c("1", "2", "3")替代useSeqnames = c(1, 2, 3)。您始終可以使用來檢查當(dāng)前R會(huì)話中是否需要染色體前綴getArchRChrPrefix()鸠匀。

1.6創(chuàng)建箭頭文件

在本教程的其余部分中蕉斜，我們將使用來自Granja *等人的造血細(xì)胞降采樣數(shù)據(jù)集的數(shù)據(jù)。自然生物技術(shù)2019。這包括來自骨髓單核細(xì)胞（BMMC）蛛勉，外周血單核細(xì)胞（PBMC）和來自骨髓的CD34 +造血干細(xì)胞和祖細(xì)胞（CD34 BMMC）的數(shù)據(jù)。

該數(shù)據(jù)以片段文件的形式下載睦柴，其中包含所有比對(duì)的測(cè)序片段的起始和終止基因組坐標(biāo)诽凌。片段文件是10x Genomics分析平臺(tái)（和其他平臺(tái)）的基本文件類型之一，可以從任何BAM文件輕松創(chuàng)建坦敌。有關(guān)創(chuàng)建自己的片段文件以輸入到ArchR的信息侣诵，請(qǐng)參見10x Genomics網(wǎng)站。

有了片段文件后狱窘，我們會(huì)將它們的路徑作為字符向量提供給createArrowFiles()杜顺。在創(chuàng)建過程中，一些基本的元數(shù)據(jù)和矩陣被添加到每個(gè)Arrow文件中蘸炸，其中包括一個(gè)“ TileMatrix”（包含跨全基因組范圍500 bp的條帶的插入計(jì)數(shù)addTileMatrix()）（請(qǐng)參見參考資料）和“ GeneScoreMatrix”躬络，該存儲(chǔ)基于在圖塊中加權(quán)插入計(jì)數(shù)的預(yù)測(cè)基因表達(dá)在基因啟動(dòng)子附近（請(qǐng)參閱參考資料addGeneScoreMatrix()）。

可以使用該getTutorialData()功能下載教程數(shù)據(jù)搭儒。教程數(shù)據(jù)的大小約為0.5 GB穷当。如果您已經(jīng)在當(dāng)前工作目錄中下載了該教程，則ArchR將繞過下載淹禾。

library(ArchR)

inputFiles <- getTutorialData("Hematopoiesis")
inputFiles

scATAC_BMMC_R1
“HemeFragments / scATAC_BMMC_R1.fragments.tsv.gz”
scATAC_CD34_BMMC_R1
“HemeFragments / scATAC_CD34_BMMC_R1.fragments.tsv.gz”
scATAC_PBMC_R1
“HemeFragments / scATAC_PBMC_R1.fragments.tsv.gz”

與往常一樣馁菜，在開始項(xiàng)目之前，我們必須為并行化設(shè)置ArchRGenome和默認(rèn)值threads铃岔。

addArchRGenome("hg19")

將默認(rèn)基因組設(shè)置為Hg19汪疮。

addArchRThreads(threads = 16) 

將并行線程的默認(rèn)數(shù)量設(shè)置為16。

現(xiàn)在毁习，我們將創(chuàng)建耗時(shí)10-15分鐘的箭頭文件智嚷。對(duì)于每個(gè)樣本，此步驟將：

1. 從提供的輸入文件中讀取可訪問的片段纺且。
2. 計(jì)算每個(gè)細(xì)胞的質(zhì)量控制信息（即TSS富集得分和核小體信息）纤勒。
3. 根據(jù)質(zhì)量控制參數(shù)過濾單元。
4. 使用500 bp的條帶創(chuàng)建全基因組的TileMatrix隆檀。
5. 使用geneAnnotation調(diào)用時(shí)定義的自定義創(chuàng)建GeneScoreMatrix addArchRGenome()摇天。

ArrowFiles <- createArrowFiles(
  inputFiles = inputFiles,
  sampleNames = names(inputFiles),
  filterTSS = 4, #Dont set this too high because you can always increase later
  filterFrags = 1000, 
  addTileMat = TRUE,
  addGeneScoreMat = TRUE
)

使用addArchRGenome（Hg19）設(shè)置的GeneAnnotation！
使用addArchRGenome（Hg19）設(shè)置的GeneAnnotation恐仑！
ArchR記錄到：ArchRLogs / ArchR-createArrows-dfa159ddbf6e-Date-2020-04-15_Time-09-21-27.log
如果有問題泉坐，請(qǐng)使用logFile向github報(bào)告！
清理臨時(shí)文件
2020-04-15 09:21:28：批量執(zhí)行w / safelapply !裳仆，已過去0分鐘腕让。
ArchR記錄成功到：ArchRLogs / ArchR-createArrows-dfa159ddbf6e-Date-2020-04-15_Time-09-21-27.log

我們可以檢查ArrowFiles對(duì)象，以查看它實(shí)際上只是Arrow文件路徑的字符向量。

ArrowFiles

“ scATAC_BMMC_R1.arrow”“ scATAC_CD34_BMMC_R1.arrow”
“ scATAC_PBMC_R1.arrow”

1.7每單元質(zhì)量控制

scATAC-seq數(shù)據(jù)的嚴(yán)格質(zhì)量控制（QC）對(duì)于消除低質(zhì)量細(xì)胞的作用至關(guān)重要纯丸。在ArchR中偏形，我們考慮數(shù)據(jù)的三個(gè)特征：

1. 獨(dú)特核片段的數(shù)量（即不映射到線粒體DNA）。
2. 信噪比觉鼻。低的信噪比通常歸因于死亡或垂死的細(xì)胞俊扭，這些細(xì)胞具有脫色的DNA，可以在全基因組范圍內(nèi)進(jìn)行隨機(jī)轉(zhuǎn)座坠陈。
3. 片段大小分布萨惑。由于核小體的周期性，我們預(yù)計(jì)會(huì)看到碎片消失仇矾，這些碎片是包裹在核小體周圍的DNA長(zhǎng)度（大約147 bp）庸蔼。

第一個(gè)度量標(biāo)準(zhǔn)是唯一的核片段，很簡(jiǎn)單-可用片段很少的細(xì)胞將無法提供足夠的數(shù)據(jù)來做出有用的解釋贮匕，因此應(yīng)排除在外姐仅。

第二個(gè)指標(biāo)，即信噪比刻盐，被計(jì)算為TSS富集得分萍嬉。傳統(tǒng)的批量ATAC-seq分析已使用此TSS富集得分作為確定信號(hào)到背景（例如，ENCODE項(xiàng)目）的標(biāo)準(zhǔn)工作流程的一部分）隙疚。我們和其他人已經(jīng)發(fā)現(xiàn)壤追，在批量ATAC-seq和scATAC-seq中測(cè)試的大多數(shù)細(xì)胞類型中，TSS富集具有代表性供屉。TSS富集得分度量標(biāo)準(zhǔn)背后的想法是行冰，由于與啟動(dòng)子結(jié)合的大型蛋白質(zhì)復(fù)合物，與其他基因組區(qū)域相比伶丐，ATAC-seq數(shù)據(jù)普遍在基因TSS區(qū)域富集悼做。通過查看以這些TSS區(qū)域?yàn)橹行牡拿總€(gè)堿基對(duì)的可及性，我們發(fā)現(xiàn)相對(duì)于側(cè)翼區(qū)域（兩個(gè)方向的遠(yuǎn)端1900-2000 bp）而言哗魂，局部富集肛走。相對(duì)于這些側(cè)翼區(qū)域，該富集峰（以TSS為中心）之間的比率表示TSS富集得分录别。

傳統(tǒng)上朽色，針對(duì)每個(gè)批量ATAC-seq樣本計(jì)算每個(gè)堿基對(duì)的可訪問性，然后使用此配置文件確定TSS富集得分组题。在scATAC-seq中以每個(gè)單元為單位執(zhí)行此操作相對(duì)較慢且計(jì)算量很大葫男。為了準(zhǔn)確估算每個(gè)單細(xì)胞的TSS富集得分，我們計(jì)算以每個(gè)單堿基TSS位置為中心的50 bp區(qū)域內(nèi)的平均可訪問性崔列，并將其除以TSS側(cè)翼位置的平均可訪問性（+/- 1900 – 2000 bp ）梢褐。該近似值與原始方法高度相關(guān)（R> 0.99），并且值在大小上非常接近。

第三個(gè)度量標(biāo)準(zhǔn)盈咳，片段大小分布通常不那么重要耿眉，但始終可以手動(dòng)檢查。由于DNA包裹核小體的模式化方式鱼响，我們希望在我們的數(shù)據(jù)中看到片段大小分布的核小體周期性鸣剪。之所以出現(xiàn)這些丘陵和山谷，是因?yàn)槠伪仨毧缭?热押、1西傀、2等斤寇。核小體（Tn5無法切割緊緊包裹在核小體周圍的DNA桶癣。

默認(rèn)情況下，在ArchR中娘锁，將通過過濾器細(xì)胞識(shí)別為那些TSS富集得分大于4且具有1000個(gè)以上獨(dú)特核片段的細(xì)胞牙寞。重要的是要注意，TSS富集得分的實(shí)際數(shù)值取決于所使用的TSS的集合莫秆。ArchR中的默認(rèn)值是為人類數(shù)據(jù)設(shè)計(jì)的间雀，在運(yùn)行時(shí)更改默認(rèn)閾值可能很重要createArrowFiles()。

創(chuàng)建Arrow文件將在當(dāng)前工作目錄中創(chuàng)建一個(gè)名為“ QualityControl”的文件夾镊屎，該文件夾將包含與每個(gè)樣本關(guān)聯(lián)的2個(gè)圖惹挟。第一張圖顯示了log10(unique nuclear fragments)vs TSS富集得分，并指出了虛線所用的閾值缝驳。第二個(gè)顯示片段大小分布连锯。

對(duì)于我們的教程數(shù)據(jù)，我們有三個(gè)示例用狱，如下所示：

對(duì)于BMMC：

對(duì)于CD34 BMMC：

參考材料：

https://www.archrproject.com/