cellranger是10X GENOMICS单细胞上游分析的软件,主要有4个流程 mkfastq、定量 count、组合 aggr、reanalyze。
如果是bcl原始测序数据,需用mkfastq转换为fastq格式(根据index将reads分配至不同的样本)。
如果是fastq格式数据,则可直接用count命令定量,得到表达矩阵,然后用aggr命令整合样本(比如实验组有多个重复样本),最后reanalyze进行后续降维聚类等等分析。(最简单的流程:如果是单个样本,只用count命令+R包即可)
本教程主要目的是从SRA或者Fastq文件完成cellranger count流程得到10x的三个文件。
流程如下
1.cellranger软件下载及安装
在网页简单注册后就可以获取wget下载地址。
首先创建并激活一个小环境
conda create -n cellranger
conda activate cellranger
下载cellranger并解压
cd /opt
wget -O cellranger-6.0.2.tar.gz "https://cf.10xgenomics.com/releases/cell-exp/cellranger-6.0.2.tar.gz?Expires=1624477084&Policy=eyJTdGF0ZW1lbnQiOlt7IlJlc291cmNlIjoiaHR0cHM6Ly9jZi4xMHhnZW5vbWljcy5jb20vcmVsZWFzZXMvY2VsbC1leHAvY2VsbHJhbmdlci02LjAuMi50YXIuZ3oiLCJDb25kaXRpb24iOnsiRGF0ZUxlc3NUaGFuIjp7IkFXUzpFcG9jaFRpbWUiOjE2MjQ0NzcwODR9fX1dfQ__&Signature=QCESOlCkbmpiDLNQEHVVCZWRUQlK01zJ28z34ezOgcb2dH7yyv0zZvWw816nE7jrOfiuiD2COZ6zf8kaL~7ndl9sfKQ~JLSWYbgZgUXb6fjehKNOJggzd32mS29lZ1cAFZBgwH~pmYEmFIIx1WIuyEKi6XZ4O6Yquc0~fUA80ZkdMoNrDGGXtgn7RgRoK4MWwGgQtsufw9J5wLXe5XQG70cmg14wd-ZGjrboK~LMBDSYfkZr2YG8Sl2ScJIbB9xKfszcyXlq65EQwFuwzmSAxvNh9uIr9YlfSeUM-uNqdc1hYkin4Q1-1nGEfAaudHgzddD45-KxBKfv0KaL-vcLBA__&Key-Pair-Id=APKAI7S6A5RYOXBWRPDA"
tar -xzvf cellranger-6.0.2.tar.gz</pre>
一般来说,软件以及配套的参考基因组都需要下载,下载速度就取决于你自己的网路情况啦,建议nohup到后台,等待即可。
添加到环境变量,方便后续使用
vim ~/.bashrc
$ export PATH=/opt/cellranger-6.0.2:$PATH
source ~/.bashrc</pre>
2.下载参考基因组并解压
wget https://cf.10xgenomics.com/supp/cell-exp/refdata-gex-GRCh38-2020-A.tar.gz
tar -xzvf refdata-gex-GRCh38-2020-A.tar.gz</pre>
3.准备数据
既然你都要学cellranger了,大概率上你已经有了SRA或者fastq数据,有了服务器,linux知识也有所了解,关于数据的下载我就不赘述了。
本次演示我们的数据来自2018年9月的NC文章Acquired cancer resistance to combination immunotherapy from transcriptional loss of class I HLA。为了展示方便,我们只使用其中一个SRR数据。
认识10x的fastq数据文件
官网给指出来了文件名规则:https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/2.0/using/fastq-input#wrongname ,如果你的fastq数据不是这样命名,就需要自行更改过来了。
zless 查看文件大小
zless -SN SRR7722937_1.fastq.gz
zless -SN SRR7722937_2.fastq.gz
zless -SN SRR7722937_3.fastq.gz</pre>
其中第一个文件的所有序列都是8bp,第二个文件都是26bp,第三个文件都是91bp,初步判断,第三个文件是测序reads
如果要理解这3个文件的区别,同理,也是需要自己去学习了解10x的原理:
看看测序时每个run cycle做了什么事:
利用illumina边合成变测序(sequencing by synthesis ,SBS),每一个cycle都是一个碱基,因此用cycle数可以表示测序长度
首先,1-26个cycle就是测序得到了26个碱基,先是16个Barcode碱基,然后是10个UMI碱基;
然后,27-34这8个cycle得到了8个碱基,就是i7的sample index;
最后35-132个cycle得到了98个碱基,就是转录本reads
4. 使用Cell Rangercount进行定量
Cell Ranger主要的流程有:拆分数据 mkfastq、细胞定量 count、定量组合 aggr、调参reanalyze,还有一些小工具比如mkref、mkgtf、upload、sitecheck、mat2csv、vdj、mkvdjref、testrun。
但是,大概率上,我们只需要使用它的定量流程,就是 cellranger count 命令。
>cellranger count --id=sample6231429 \
--transcriptome=/home/rstudio/data/ref/refdata-cellranger-GRCh38-1.2.0 \
--fastqs=/home/rstudio/data/raw0 \
--sample=SRR7722937 \
--nosecondary
# id指定输出文件存放目录名
# transcriptome指定与CellRanger兼容的参考基因组
# fastqs指定mkfastq或者自定义的测序文件
# sample要和fastq文件的前缀中的sample保持一致,作为软件识别的标志
# expect-cells指定复现的细胞数量,这个要和实验设计结合起来
# nosecondary 只获得表达矩阵,不进行后续的降维、聚类和可视化分析(因为后期会自行用R包去做)
服务器配置不一样,这个cellranger count流程运行时间不一样,这一个样本是约3G的fq文件数据走这个流程是2小时。
我的服务器是配置是32核40线程128G内存
当然你也可以写个脚本转后台,这都是后话了。
5.输出结果
运行结束得到结果如下:
Outputs:
- Run summary HTML: /home/rstudio/data/sample6231429/outs/web_summary.html
- Run summary CSV: /home/rstudio/data/sample6231429/outs/metrics_summary.csv
- BAM: /home/rstudio/data/sample6231429/outs/possorted_genome_bam.bam
- BAM index: /home/rstudio/data/sample6231429/outs/possorted_genome_bam.bam.bai
- Filtered gene-barcode matrices MEX: /home/rstudio/data/sample6231429/outs/filtered_gene_bc_matrices
- Filtered gene-barcode matrices HDF5: /home/rstudio/data/sample6231429/outs/filtered_gene_bc_matrices_h5.h5
- Unfiltered gene-barcode matrices MEX: /home/rstudio/data/sample6231429/outs/raw_gene_bc_matrices
- Unfiltered gene-barcode matrices HDF5: /home/rstudio/data/sample6231429/outs/raw_gene_bc_matrices_h5.h5
- Secondary analysis output CSV: null
- Per-molecule read information: /home/rstudio/data/sample6231429/outs/molecule_info.h5
- Loupe Cell Browser file: null
Pipestance completed successfully!
Saving pipestance info to sample6231429/sample6231429.mri.tgz
(cellrangerze) root 08:42:19 /home/rstudio/data
输出结果说明:
filtered_gene_bc_matrices:是重要的一个目录,下面又包含了 barcodes.tsv、features.tsv、matrix.mtx,是下游Seurat、Scater、Monocle等分析的输入文件。
web_summary.html:质控比对报告(一般认为外显子的比对率要在60%以上)。barcode用来标记细胞,UMI用来标记转录本;其次,barcodes数量时要大于细胞数量的(以保证每个细胞都会有barcode来进行区分)。
6.一些术语解释
index标记样本,barcode标记细胞,UMI标记转录本。
i7 sample index (library barcode)
是加到Illumina测序接头上的,保证多个测序文库可以在同一个flow-cell上或者同一个lane上进行混合测序(multiplexed)。它的作用就是在CellRanger的mkfastq 功能中体现出来的,它自动识别样本index名称(例如:SA-GA-A1),将具有相同4种oligo的fq文件组合在一起,表示同一个样本。它保证了一个测序lane上可以容纳多个样本。一个index set有4个oligos。
Barcode
是10X特有的,用来区分GEMs,也就是对细胞做了一个标记。一般在拆分混样测序数据(demultiplexing)这个过程后进行操作,当然这也很符合原文的操作。
UMI
UMI就是Unique Molecular Identifier,由4-10个随机核苷酸组成,在mRNA反转录后,进入到文库中,每一个mRNA随机连上一个UMI,根据PCR结果可以计数不同的UMI,最终统计mRNA的数量。它的主要作用是,处理PCR 扩增偏差,因为起始文库很小时需要的PCR扩增次数就越多,因为越容易引入扩增误差。
v2 chemistry和v3 chemistry的区别
和V2相比,V3试剂盒中所用的UMI和PolyT的长度都发生了变化,从而导致测序得到的R1和R2端的序列长度也不一致,V2试剂盒的R1端长度为26bp, 包含16bp的barcode和10bp的UMI序列,V3试剂盒的R1端长度为28bp, 包含16bp的barcode和12bp的UMI序列;V2试剂盒的R2端为98bp, V3试剂盒的R2端为91bp。
参考
cellranger更新到4啦(全新使用教程) | 生信菜鸟团 (bio-info-trainee.com)
cellranger更新到5啦(全新使用教程) - 云+社区 - 腾讯云 (tencent.com)
CellRanger使用学习 - 简书 (jianshu.com)