在分析单细胞转录组测序数据前,我们需要将原始的fastq测序文件转化成我们可读的数据,最好的方法是将数据整理为基因-细胞表达矩阵,如此一来后续的分析就能够直接读取到某细胞的某基因表达量是多少。
梗概
本文为整理fastq文件到基因表达水平的方法,为下游分析做准备。
软件安装与数据下载
- Cell Ranger软件(需要linux系统)
下载地址:https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/installation
如何使用:https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/what-is-cell-ranger - 基因组参考数据
下载地址:https://support.10xgenomics.com/single-cell-gene-expression/software/downloads/latest
在将测序数据的序列匹对到一个特定的基因时,我们需要一个该物种的全基因组序列参考数据,起到一个“字典”的作用,告诉软件怎么样的序列是属于什么基因。人类的数据一般采用GRCh38或hg19。 - 分析数据
这里采用的是Human Cell Atlas计划中的脾脏组织测序数据。
下载地址:
https://preview.data.humancellatlas.org
Human Cell Atlas Preview Datasets数据为Human Cell Atlas(HCA)的第一批公开的单细胞测序数据库,HCA计划目的在于收集人类各个器官组织的单细胞测序数据,汇合成一个完整的人类单细胞图谱,便于后续人们对单个的组织的单细胞测序数据研究时的横向对比和参照。
目前一共有三个数据,分别为人类免疫细胞的普查数据、缺血敏感性人类脾脏组织和黑色素瘤浸润机制基质和免疫细胞的单细胞RNA-seq数据。
其中前两个数据集是通过10x genomics公司的测序方法测得,为fastq格式,正好可用于cell ranger进行整理和分析。
数据整理
在完成上述的准备之后,整理过程十分简单,只需要一个命令。运行这一个命令就可以将fastq文件解析为基因-细胞表达矩阵,并且可以得到一些初步的细胞表达统计结果。结果包括:
- 细胞数量、基因匹对率等基础信息
- 降维分析,将细胞投射到二维空间(t-SNE)
- 自动聚类分析,将具有相似表达谱的细胞组合在一起
- 在所选cluster之间差异表达的基因列表
- 显示测序深度减少对观察到的文库复杂性的影响
- 显示测序深度减少对检测到的中值基因的影响
## 利用cd命令进入cell ranger软件所在的文件夹
## cellranger count命令
## --id给你这次的运行七个名字,如sample345
## --fastqs 输入分析数据所在路径
## --transcriptome 输入参考基因组所在路径
cellranger count --id=HCATisStabAug177078016 \
> --fastqs=/ischaemic_sensitivity/17a3d288-01a0-464a-9599-7375fda3353d/ \
> --transcriptome=/cellranger/
……
2018-09-20 19:00:41 [runtime] (update) ID.HCATisStabAug177078016.SC_RNA_COUNTER_CS.SC_RNA_COUNTER._BASIC_SC_RNA_COUNTER.EXTRACT_READS.fork0 chunks running (0/41 completed)
2018-09-20 19:06:33 [runtime] (update) ID.HCATisStabAug177078016.SC_RNA_COUNTER_CS.SC_RNA_COUNTER._BASIC_SC_RNA_COUNTER.EXTRACT_READS.fork0 chunks running (1/41 completed)
2018-09-20 19:12:42 [runtime] (update) ID.HCATisStabAug177078016.SC_RNA_COUNTER_CS.SC_RNA_COUNTER._BASIC_SC_RNA_COUNTER.EXTRACT_READS.fork0 chunks running (2/41 completed)
2018-09-20 19:18:40 [runtime] (update) ID.HCATisStabAug177078016.SC_RNA_COUNTER_CS.SC_RNA_COUNTER._BASIC_SC_RNA_COUNTER.EXTRACT_READS.fork0 chunks running (3/41 completed)
2018-09-20 19:24:40 [runtime] (update) ID.HCATisStabAug177078016.SC_RNA_COUNTER_CS.SC_RNA_COUNTER._BASIC_SC_RNA_COUNTER.EXTRACT_READS.fork0 chunks running (4/41 completed)
2018-09-20 19:30:41 [runtime] (update) ID.HCATisStabAug177078016.SC_RNA_COUNTER_CS.SC_RNA_COUNTER._BASIC_SC_RNA_COUNTER.EXTRACT_READS.fork0 chunks running (4/41 completed)
需要耐心等待,这个命令运行完需要5个小时。
运行完成后,会给出输出文件:
Outputs:
- Run summary HTML: /outs/web_summary.html
- Run summary CSV: /outs/metrics_summary.csv
- BAM: /outs/possorted_genome_bam.bam
- BAM index: /outs/possorted_genome_bam.bam.bai
- Filtered gene-barcode matrices MEX: /outs/filtered_gene_bc_matrices
- Filtered gene-barcode matrices HDF5: /outs/filtered_gene_bc_matrices_h5.h5
- Unfiltered gene-barcode matrices MEX: /outs/raw_gene_bc_matrices
- Unfiltered gene-barcode matrices HDF5: /outs/raw_gene_bc_matrices_h5.h5
- Secondary analysis output CSV: /outs/analysis
- Per-molecule read information: /outs/molecule_info.h5
- Loupe Cell Browser file: /outs/cloupe.cloupe
Waiting 6 seconds for UI to do final refresh.
Pipestance completed successfully!
Saving pipestance info to HCATisStabAug177078016/HCATisStabAug177078016.mri.tgz
输出在cellranger文件夹中产生了以sampleID命名的文件夹,也就是HCATisStabAug177078016。
将web_summary.html下载到本地查看
可以通过单击左上角的“Summary”查看运行摘要。描述了检测细胞的测序质量和各种特征。
这里展示的是匹对率,如果匹对率太低,有可能是测序质量不佳,也有可能是选错参考基因组。
点击左上角的“Analysis”即可查看自动二次分析结果。二次分析提供以下内容:
- 降维分析,将细胞投射到二维空间(t-SNE)
- 自动聚类分析,将具有相似表达谱的细胞组合在一起
- 在所选cluster之间差异表达的基因列表
- 显示测序深度减少对观察到的文库复杂性的影响
- 显示测序深度减少对检测到的中值基因的影响
这里显示的是每个细胞条形码的总UMI计数。每个点表示一个细胞,颜色表示UMI含量。具有较大UMI计数的细胞可能具有比具有较少UMI计数的细胞更高的RNA含量,也就是越红的细胞RNA含量越高。坐标轴对应于由t-SNE算法产生的二维嵌入。在该空间中,彼此接近的细胞对具有比彼此远离的细胞更相似的基因表达谱。显示器限于10000个单元的随机子集。
这是通过自动聚类算法对每个细胞条形码进行聚类。聚类将具有相似表达谱的单元组合在一起。坐标轴对应于由t-SNE算法产生的二维嵌入。在该空间中,彼此接近的细胞对具有比彼此远离的细胞更相似的基因表达谱。显示器限于10000个单元的随机子集,K-均值最大为K = 20。可使用Loupe(tm)Cell Browser查看整个数据集。