引文

上期推文【scATAC-seq3:常用工具—SnapATAC簡介】當(dāng)中，我們主要對SnapATAC這一個(gè)工具的特點(diǎn)進(jìn)行了簡單的介紹。在本期推文當(dāng)中，我們將繼續(xù)上一次的話題，簡單介紹scATAC-seq的上游分析流程，即最常用的Cellranger和用于SnapATAC分析的上游分析軟件snaptools。

Cellranger 上游分析

1）版本的選擇

對于Cellranger ATAC的版本相比于RNA而言要少很多，主要可以分為2.0和1.2及之前的版本。2.0版本相比于1.2之前的版本，在算法方面有了比較大的改動(dòng)。

首先針對于標(biāo)記PCR重復(fù)這一流程，1.2之前的版本主要以起始位置和末端位置為基礎(chǔ)進(jìn)行標(biāo)記，造成的結(jié)果是序列的重復(fù)率會隨著可及性的增加而增加。2.0版本則是除了基于起始和末端位置以外，同時(shí)根據(jù)散列的barcode進(jìn)行標(biāo)記，能夠提高對標(biāo)記重復(fù)的準(zhǔn)確度。

此外，新舊版本的差異主要體現(xiàn)在peak calling。在舊版本當(dāng)中，peak calling主要是基于計(jì)算得到的全局閾值，即全局閾值以上的含有平滑信號的連續(xù)區(qū)域，因此并不能準(zhǔn)確識別所有的motif位點(diǎn)。新版本中對背景噪聲更加敏感，準(zhǔn)確度更高。

2）建立索引

Cellranger ATAC的建立索引主要需要三個(gè)文件：

• 參考基因組文件、
• GENCODE上的功能元件注釋文件、
• 轉(zhuǎn)錄因子及其motif文件。

以建立人的GRCh38的索引為例，則需要:

• Homo_sapiens.GRCh38.dna.primary_assembly.fa.gz、
• gencode.v32.primary_assembly.annotation.gtf.gz、
• JASPAR2018_CORE_non-redundant_pfms_jaspar.txt

這三個(gè)文件為基礎(chǔ)進(jìn)行建立。

具體建立索引的步驟可以參考 https://support./single-cell-atac/software/release-notes/references

3）cellranger-atac count

scATAC的現(xiàn)在基本上從公司拿到的數(shù)據(jù)都是fastq結(jié)尾的原始文件，則直接可以從cellranger-atac count這個(gè)步驟開始運(yùn)行。

cellranger-atac count   --id=sample345 \
                        --reference=/opt/refdata-cellranger-arc-GRCh38-2020-A-2.0.0 \
                        --fastqs=/home/jdoe/runs/HAWT7ADXX/outs/fastq_path \
                        --sample=mysample \
                        --localcores=8 \
                        --localmem=64

和RNA相似，只是將參數(shù)中的--transcriptome調(diào)整為--reference。需要注意的是如果沒有設(shè)置localcores和localmem，將會運(yùn)用系統(tǒng)中可用的所有線程和內(nèi)存。

4) 運(yùn)行結(jié)果

Outputs:
- Per-barcode fragment counts & metrics:        /home/jdoe/runs/sample345/outs/singlecell.csv
- Position sorted BAM file:                     /home/jdoe/runs/sample345/outs/possorted_bam.bam
- Position sorted BAM index:                    /home/jdoe/runs/sample345/outs/possorted_bam.bam.bai
- Summary of all data metrics:                  /home/jdoe/runs/sample345/outs/summary.json
- HTML file summarizing data & analysis:        /home/jdoe/runs/sample345/outs/web_summary.html
- Bed file of all called peak locations:        /home/jdoe/runs/sample345/outs/peaks.bed
- Raw peak barcode matrix in hdf5 format:       /home/jdoe/runs/sample345/outs/raw_peak_bc_matrix.h5
- Raw peak barcode matrix in mex format:        /home/jdoe/runs/sample345/outs/raw_peak_bc_matrix
- Directory of analysis files:                  /home/jdoe/runs/sample345/outs/analysis
- Filtered peak barcode matrix in hdf5 format:  /home/jdoe/runs/sample345/outs/filtered_peak_bc_matrix.h5
- Filtered peak barcode matrix in mex format:   /home/jdoe/runs/sample345/outs/filtered_peak_bc_matrix
- Barcoded and aligned fragment file:           /home/jdoe/runs/sample345/outs/fragments.tsv.gz
- Fragment file index:                          /home/jdoe/runs/sample345/outs/fragments.tsv.gz.tbi
- Filtered tf barcode matrix in hdf5 format:    /home/jdoe/runs/sample345/outs/filtered_tf_bc_matrix.h5
- Filtered tf barcode matrix in mex format:     /home/jdoe/runs/sample345/outs/filtered_tf_bc_matrix
- Loupe Browser input file:                     /home/jdoe/runs/sample345/outs/cloupe.cloupe
- csv summarizing important metrics and values: /home/jdoe/runs/sample345/outs/summary.csv
- Annotation of peaks with genes:               /home/jdoe/runs/sample345/outs/peak_annotation.tsv
- Peak-motif associations:                      /home/jdoe/runs/sample345/outs/peak_motif_mapping.bed

對于不同的下游分析軟件，讀取的文件是不同的。

• 是ArchR，讀取的是fragments.tsv.gz文件；
• 是SnapATAC，推薦的方式是通過將bam文件進(jìn)行轉(zhuǎn)化為snap文件或者也可以通過fragments.tsv.gz文件產(chǎn)生snap文件；
• Signac則是需要singlecell.csv、filtered_peak_bc_matrix.h5、fragments.tsv.gz三個(gè)文件為基礎(chǔ)進(jìn)行讀取。

所以，我們經(jīng)常出現(xiàn)的情況是ArchR讀取的細(xì)胞數(shù)量和Cellranger產(chǎn)生的summary中的細(xì)胞數(shù)量是不同的。

snaptools上游分析

上游分析流程（建立在fastq基礎(chǔ)上）主要含有五個(gè)步驟：

1）測序文庫拆分

2）建立索引文件

3）比對

4）數(shù)據(jù)預(yù)處理

5）產(chǎn)生表達(dá)矩陣

對于第一步主要是通過python進(jìn)行實(shí)現(xiàn)，可以參考作者提供的代碼 https://github.com/r3fang/SnapTools/blob/master/snaptools/dex_fastq.py 。其余的步驟（2-5）可以參考https://github.com/r3fang/SnapATAC/wiki/FAQs。如果之前運(yùn)行過Cellranger，則可以通過產(chǎn)生的bam文件進(jìn)行轉(zhuǎn)換。

總結(jié)

本期我們主要是簡單介紹了一下Cellranger ATAC的上游分析流程。總的來說，Cellranger ATAC的運(yùn)行時(shí)間相比RNA運(yùn)行的時(shí)間更長，而在下游分析的過程當(dāng)中也發(fā)現(xiàn)scATAC-seq相比于scRNA-seq的運(yùn)行時(shí)間和內(nèi)存需要的更多。在下一期推文當(dāng)中，我們會開始介紹scATAC-seq的下游分析流程。