公共數(shù)據(jù)庫中的SRA 單細胞轉錄組數(shù)據(jù)究竟包含了哪些數(shù)據(jù)？CellRanger怎樣利用10x平臺下機數(shù)據(jù)進行下游一系列分析？這篇文章簡單記錄CellRanger 包括的主要分析步驟，純理論。
SRA 原始數(shù)據(jù)轉fastq
公共數(shù)據(jù)庫的SRA 數(shù)據(jù)需要借助fastq-dump 轉為fastq文件集索，然后進行質(zhì)控、CellRanger定量等操作汇跨。相較于普通轉錄組數(shù)據(jù)务荆，原始SRA數(shù)據(jù)會得到3個fastq文件，分別是Library 的Index(8bp)文件穷遂，包括長度為26bp 的Barcode（16bp）和UMI(10bp)的Read1文件函匕，和測序reads文件。
conda install -c bioconda sra-tools ## 安裝軟件
wkd=/home/project/single-cell/MCC
cd $wkd/raw/P2586-4
cat SRR_Acc_List-2586-4.txt |while read i
do
time fastq-dump --gzip --split-files -A $i ${i}.sra && echo "** ${i}.sra to fastq done **"
done
### 單細胞數(shù)據(jù)參數(shù)為 --split-files 而不是 --split-3
i7 sample index (library barcode）
是加到Illumina測序接頭上的蚪黑，保證多個測序文庫可以在同一個flow-cell上或者同一個lane上進行混合測序(multiplexed)盅惜。它的作用就是在CellRanger的mkfastq 功能中體現(xiàn)出來的，它自動識別樣本index名稱（例如：SA-GA-A1）忌穿，將具有相同4種oligo的fq文件組合在一起抒寂，表示同一個樣本。它保證了一個測序lane上可以容納多個樣本
Barcode
是10X特有的掠剑，用來區(qū)分GEMs屈芜，也就是對細胞做了一個標記。一般在拆分混樣測序數(shù)據(jù)（demultiplexing）這個過程后進行操作，當然這也很符合原文的操作井佑。
UMI
UMI就是Unique Molecular Identifier属铁，由4-10個隨機核苷酸組成，在mRNA反轉錄后躬翁，進入到文庫中焦蘑，每一個mRNA隨機連上一個UMI，根據(jù)PCR結果可以計數(shù)不同的UMI姆另，最終統(tǒng)計mRNA的數(shù)量喇肋。它的主要作用是，處理PCR 擴增偏差迹辐，因為起始文庫很小時需要的PCR擴增次數(shù)就越多蝶防，因為越容易引入擴增誤差。
fastq文件更名
為什么要更名明吩？CellRanger 定量過程輸入文件指定命名格式间学。
如何更名？下圖格式：
# 比如印荔，將原來的SRR7692286_1.fastq.gz改成SRR7692286_S1_L001_I1_001.fastq.gz
# 依次類推低葫，將原來_2的改成R1，將_3改成R2
cat SRR_Acc_List-9245-3.txt | while read i ;do (mv ${i}_1*.gz ${i}_S1_L001_I1_001.fastq.gz;mv ${i}_2*.gz ${i}_S1_L001_R1_001.fastq.gz;mv ${i}_3*.gz ${i}_S1_L001_R2_001.fastq.gz);done
fastq 文件質(zhì)控
主要查看Q30比例仍律。
# 以P2586-4為例
mkdir -p $wkd/qc
cd $wkd/qc
find $wkd/raw/P2586-4 -name '*R1*.gz'>P2586-4-id-1.txt
find $wkd/raw/P2586-4 -name '*R2*.gz'>P2586-4-id-2.txt
cat P2586-4-id-1.txt P2586-4-id-2.txt >P2586-4-id-all.txt
cat P2586-4-id-all.txt| xargs fastqc -t 20 -o ./
CellRanger 流程
它主要包括四個主要基因表達分析流程：
-
cellranger mkfastq ： 它借鑒了Illumina的
bcl2fastq嘿悬，可以將一個或多個lane中的混樣測序樣本按照index標簽生成樣本對應的fastq文件。即對下機數(shù)據(jù)base calling files轉為fastq文件水泉。 -
cellranger count ：利用
mkfastq生成的fq文件善涨，進行比對(基于STAR)、過濾草则、UMI計數(shù)钢拧。利用細胞的barcode生成gene-barcode矩陣，然后進行樣本分群炕横、基因表達分析 -
cellranger aggr ：接受
cellranger count的輸出數(shù)據(jù)源内，將同一組的不同測序樣本的表達矩陣整合在一起，比如tumor組原來有4個樣本份殿，PBMC組有兩個樣本膜钓，現(xiàn)在可以使用aggr生成最后的tumor和PBMC兩個矩陣，并且進行標準化去掉測序深度的影響 -
cellranger reanalyze ：接受
cellranger count或cellranger aggr生成的gene-barcode矩陣卿嘲，使用不同的參數(shù)進行降維颂斜、聚類。屬于定制化分析腔寡。
它的結果主要是包含有細胞信息的BAM, MEX, CSV, HDF5 and HTML文件
CellRanger 軟件安裝及測試實戰(zhàn)操作參考：CellRanger走起(四) CellRanger流程概覽焚鲜。包括安裝詳細指南及其他一些小技巧掌唾。
拆分數(shù)據(jù) mkfastq放前、細胞定量 count忿磅、定量組合 aggr
mkfastq 拆分
目的：將每個flowcell 的Illumina sequencer's base call files (BCLs)轉為fastq文件
特色： 它借鑒了Illumina出品的bcl2fastq，另外增加了：
- 將10X 樣本index名稱與四種寡核苷酸對應起來凭语，比如A1孔是樣本SI-GA-A1葱她，然后對應的寡核苷酸是GGTTTACT, CTAAACGG, TCGGCGTC, and AACCGTAA ，那么程序就會去index文件中將存在這四種寡核苷酸的fastq組合到A1這個樣本
- 提供質(zhì)控結果似扔，包括barcode 質(zhì)量吨些、總體測序質(zhì)量如Q30、R1和R2的Q30堿基占比炒辉、測序reads數(shù)等
- 可以使用10X簡化版的樣本信息表
### 兩種使用方式
# 第一種
$ cellranger mkfastq --id=bcl \
--run=/path/to/bcl \
--samplesheet=samplesheet-1.2.0.csv
# 第二種
$ cellranger mkfastq --id=bcl \
--run=/path/to/bcl \
--csv=simple-1.2.0.csv
# 其中id指定輸出目錄的名稱豪墅，run指的是下機的原始BCL文件目錄
# 重要的就是測序lane、樣本名稱黔寇、index等信息
參考文章CellRanger走起(四) CellRanger流程概覽 中解釋了幾種不同方式輸出fastq文件的不同存放目錄偶器，可根據(jù)實際分析自行查找。
count 定量
這個過程是最重要的缝裤，它完成細胞與基因的定量屏轰，它將比對、質(zhì)控憋飞、定量都包裝了起來霎苗，內(nèi)部流程很多，但使用很簡單榛做。每個版本要求的參數(shù)是不同的唁盏，尤其是V2與V3版本存在較大差異，這里先對V2進行了解瘤睹。
# 這是示例升敲，不是真實數(shù)據(jù) #
cellranger count --id=sample345 \
--transcriptome=/opt/refdata-cellranger-GRCh38-1.2.0 \
--fastqs=/home/scRNA/runs/HAWT7ADXX/outs/fastq_path \
--sample=mysample \
--expect-cells=1000 \
--nosecondary
# id指定輸出文件存放目錄名
# transcriptome指定與CellRanger兼容的參考基因組
# fastqs指定mkfastq或者自定義的測序文件
# sample要和fastq文件的前綴中的sample保持一致，作為軟件識別的標志
# expect-cells指定復現(xiàn)的細胞數(shù)量轰传，這個要和實驗設計結合起來
# nosecondary 只獲得表達矩陣驴党，不進行后續(xù)的降維、聚類和可視化分析(因為后期會自行用R包去做)
###輸出文件
Outputs:
- Run summary HTML: /opt/sample345/outs/web_summary.html
- Run summary CSV: /opt/sample345/outs/metrics_summary.csv
- BAM: /opt/sample345/outs/possorted_genome_bam.bam
- BAM index: /opt/sample345/outs/possorted_genome_bam.bam.bai
- Filtered gene-barcode matrices MEX: /opt/sample345/outs/filtered_gene_bc_matrices
- Filtered gene-barcode matrices HDF5: /opt/sample345/outs/filtered_gene_bc_matrices_h5.h5
- Unfiltered gene-barcode matrices MEX: /opt/sample345/outs/raw_gene_bc_matrices
- Unfiltered gene-barcode matrices HDF5: /opt/sample345/outs/raw_gene_bc_matrices_h5.h5
- Secondary analysis output CSV: /opt/sample345/outs/analysis
- Per-molecule read information: /opt/sample345/outs/molecule_info.h5
- Loupe Cell Browser file: /opt/sample345/outs/cloupe.cloupe
Pipestance completed successfully!
輸出文件從上到下依次來看：
- web_summary.html：官方說明 summary HTML file
- metrics_summary.csv：CSV格式數(shù)據(jù)摘要
- possorted_genome_bam.bam：比對文件
- possorted_genome_bam.bam.bai：索引文件
- filtered_gene_bc_matrices：是重要的一個目錄获茬，下面又包含了 barcodes.tsv.gz港庄、features.tsv.gz、matrix.mtx.gz恕曲，是下游Seurat鹏氧、Scater、Monocle等分析的輸入文件
- filtered_feature_bc_matrix.h5：過濾掉的barcode信息HDF5 format
- raw_feature_bc_matrix：原始barcode信息
- raw_feature_bc_matrix.h5：原始barcode信息HDF5 format
- analysis：數(shù)據(jù)分析目錄佩谣，下面又包含聚類clustering（有graph-based & k-means）把还、差異分析diffexp、主成分線性降維分析pca、非線性降維tsne
- molecule_info.h5：下面進行aggregate使用的文件
- cloupe.cloupe：官方可視化工具Loupe Cell Browser 輸入文件
count 定量的重點和難點在于參考基因組索引的構建吊履，可詳細參考CellRanger走起(四) CellRanger流程概覽安皱，文章中還介紹了怎樣自行構建參考基因組索引，及一比對過程的重要細節(jié)信息艇炎。
aggr
當處理多個生物學樣本或者一個樣本存在多個重復/文庫時酌伊，最好的操作就是先分別對每個文庫進行單獨的count定量，然后將定量結果利用aggr組合起來.
- 得到count結果
- 構建Aggregation CSV
- 運行aggr
### step1
$ cellranger count --id=LV123 ...
... wait for pipeline to finish ...
$ cellranger count --id=LB456 ...
... wait for pipeline to finish ...
$ cellranger count --id=LP789 ...
... wait for pipeline to finish ...
## step2
# AGG123_libraries.csv
library_id,molecule_h5
LV123,/opt/runs/LV123/outs/molecule_info.h5
LB456,/opt/runs/LB456/outs/molecule_info.h5
LP789,/opt/runs/LP789/outs/molecule_info.h5
# 其中
# molecule_h5：文件molecule_info.h5 file的路徑
### step3
cellranger aggr --id=AGG123 \
--csv=AGG123_libraries.csv \
--normalize=mapped
# 結果輸出到AGG123這個目錄中
cellranger count的結果
count結果一般放在out目錄下缀踪，主要有summary和analysis兩大類
Summary
是一個HTML格式文件居砖，包括實驗捕獲的細胞數(shù)目、檢測到的基因數(shù)目驴娃、測序數(shù)據(jù)的產(chǎn)出與質(zhì)量統(tǒng)計奏候、參考基因組的比對情況。
幾個指標可以關注一下：
左上部分中唇敞，包括了reads數(shù)鼻由、barcodes數(shù)、UMI厚棵、index蕉世、Q30等統(tǒng)計值
左下是reads比對的比例，包括基因間區(qū)婆硬、外顯子狠轻、內(nèi)含子，如果比對率太低(一般認為外顯子的比對率要在60%以上)
右上圖是利用barcodes上的UMI標簽分布來估計細胞數(shù)彬犯，綠/藍色表示細胞向楼，灰色表示背景，其中Y軸表示每個barcode對應的UMI數(shù)量谐区，X軸是一定數(shù)量的UMI序列所對應barcode的數(shù)量湖蜕，比如上圖中有1000個barcodes含有10k個UMI，細胞過濾就是通過這個圖來展示的宋列。
首先明確昭抒，barcode用來區(qū)分細胞，UMI用來區(qū)分轉錄本炼杖；其次灭返，barcodes數(shù)量時要大于細胞數(shù)量的（以保證每個細胞都會有barcode來進行區(qū)分）如果圖線陡降說明
另參考文章CellRanger走起(二) CellRanger out 結果 中還記錄了一個手動構建物種GTF文件的案例嘗試，及相應的檢驗方法坤邪，可以根據(jù)需要自行學習熙含。
總之，無論是從公共數(shù)據(jù)庫下載單細胞SRA數(shù)據(jù)艇纺，還是來自10x 平臺的原始下機數(shù)據(jù)怎静，或者先經(jīng)過fastq-dump轉為fastq文件邮弹，然后質(zhì)控，然后CellRanger count 定量蚓聘，或者首先自行拆分數(shù)據(jù)mkfastq肠鲫，然后構建索引進行細胞定量，及可選的定量組合等或粮。CellRanger的總體流程見本文記錄，但是其中的各個細節(jié)仍需要自行實踐總結捞高。
參考1：CellRanger走起(二) 使用前注意事項
參考2：CellRanger走起(三) 使用初探
參考3：CellRanger走起(四) CellRanger流程概覽
參考4：CellRanger走起(二) CellRanger out 結果
