10x 空间转录组必备分析软件--Space Ranger

2022-06-30 19:47:32, 星标关注上海欧易生物医学科技有限公司

Space Ranger 软件是 10x Genomics 官方提供的空间转录组分析的配套分析软件，是一套利用亮场和荧光显微镜图像处理 Visium 空间基因表达数据的分析流程，也是在 10x Visium 空间转录组的分析中重要的分析软件。下面小欧将从以下几个方面来详细介绍：

什么是 Space Ranger

Space Ranger 包括 5 个与 Visium 空间基因表达实验相关的流程。

1. spaceranger mkfastq：其功能为将 Illumina 测序仪产生的 raw base call (BCL) 文件解析成 FASTQ 文件。一般我们拿到准备分析的数据已经是 FASTQ 序列文件，因此，不需要再运行 spaceranger mkfastq。

2. spaceranger count：其功能是从 spaceranger mkfastq 中获取 FASTQ 文件和显微镜载玻片图像，并进行对齐、组织检测、基准检测和条形码/UMI 计数。该流程使用 Visium 空间条形码生成特征点矩阵，确定聚类，并执行基因表达分析。

3. spaceranger aggr：其功能为将多个样本的 spaceranger count 产生的数据进行整合、标准化，并可以对整合后的数据进行分析。

4. spaceranger targeted-compare：其功能为将起始输入库(称为父库)与其相应的 Targeted Gene Expression 数据集进行比较。与只知道目标数据的情况相比，spaceranger targeted-compare 可以更准确地评估目标性能。提供了质量控制指标，以验证目标基因的富集程度和亲本数据的恢复。该流程仅支持新鲜冷冻组织。

5. spaceranger targeted-depth：在假设的靶向基因表达实验的背景下，总结了整个转录组分析(WTA)数据集。给定现有的 WTA 数据集和目标 panel 的 CSV 文件，spaceranger targeted-depth 计算 panel 中映射到目标基因的读取部分。该流程仅支持新鲜冷冻组织。

10x 空间转录组技术原理

由于 Space Ranger 分析与技术原理息息相关，不懂技术原理，分析起来就是一头雾水。接下来，我们一起来看看其技术原理是什么样的，如何获取原位信息呢？

10x Genomics 空间转录组用于文库构建的每张载玻片上有四个捕获区域，每个捕获区域的大小为 6.5x6.5 mm，包含近 5000 个被条形码标记的点（barcoded Spots），每个点的直径为 55 μm，点和点之间中心的距离为 100 μm，并且每个点都有一个唯一的 barcode 序列。组织切片的细胞中会释放出 mRNA，迁移到每个 Spot 的 mRNA 会被标记上相应的 barcode 序列，然后进行文库构建并进行测序。最后，根据数据的条形码信息对数据进行分析，以确定哪些数据来自哪个位置，从而实现空间基因表达的可视化。

接下来我们来进入正题，来看看 Space Ranger 是如何分析空间基因表达的。

Space Ranger 的核心原理

比起 Cell Ranger 的分析原理，Space Ranger 新增了利用成像算法来识别基因表达的具体位置步骤。

成像算法

Space Ranger 依靠图像处理算法来解决与玻片（slide ）图像相关的两个关键问题：

1. 将矩阵 barcode 和切片照片相关联

• 基准点对齐，定位每个 spot 在图像中的位置。

• 自动检测失败可以使用 Loupe Browser 手动调整。

由于每个用户在 Visium 捕获区域成像时可能会因为视角不同而设置不同，数据需要基准对齐，这样 Space Ranger 就可以知道图像中单个条形码 spot 的位置。

2. 区分组织 spot 和背景 spot

• 检测将用于下游分析的组织位置。

• 使用光学显微镜的明场图像，标记高亮度的像素点为背景，低亮度的像素点为组织。

序列比对和基因定量，生成矩阵信息

这部分和 Cell Ranger 检测基因表达原理类似，此处不展开阐述，感兴趣参见Cell Ranger 知多少？（上）

spaceranger count 使用方法

Space Ranger 软件最重要的 pipeline 是 spaceranger count。spaceranger count 运行时图片对齐方式有两种，一种软件自动识别图片进行对齐，另外一种就是先用 Loupe Browser 软件手动对齐，生成对应 json 文件提供给软件。

自动对齐：

spaceranger count --id=Sample_A \\
                   --transcriptome=spaceranger-refdata/refdata-gex-GRCh38-2020-A \\
                   --fastqs=raw_data/Sample_A/Sample_A_fastqs \\
                   --sample=Sample_A \\
                   --image=raw_data/images/Sample_A.jpg \\
                   --slide=V19J01-123 \\
                   --area=A1 \\
                   --localcores=8 \\
                   --localmem=64

手动对齐:

spaceranger count --id=Sample_A \\
                   --transcriptome=spaceranger-refdata/refdata-gex-GRCh38-2020-A \\
                   --fastqs=raw_data/Sample_A/Sample_A_fastqs \\
                   --sample=Sample_A \\
                   --image=raw_data/images/Sample_A.jpg \\
                   --slide=V19J01-123 \\
                   --area=A1 \\
                   --loupe-alignment=Sample_A.json \\
                   --localcores=8 \\
                   --localmem=64

重要参数说明：

--id：结果输出文件夹名称

--transcriptome：基因组目录

--fastqs：fastq 文件目录

--sample：原始样本名

--image：组织 H&E 染色图片，可以是 jpg 或者 tiff 格式。

--slide：使用的 10x 芯片型号

--area：样本所在芯片的区域（四通道芯片位置从上到下分别为 A1、B1、C1、D1）

--loupe-alignmen：图片手动对齐生成的 json 文件，如果是自动对齐不要需要此参数

--localcores: 限制 spaceranger 一次最多使用的核数

--localmem: 限制 spaceranger 使用的最大内存量

--slidefile: 如果要在无法访问 Internet 的环境中运行，流程将需要通过 --slidefile 参数提供 Visium 玻片布局文件。如 slide 序列号为：V19J01-123, 则在无法访问网络的环境下需要提前提供 下载好 V10L21-109.gpr。

spaceranger count 结果解读

一个成功的结果输出是这样的：

Outputs:
- Run summary HTML:                         Sample_A/outs/web_summary.html
- Outputs of spatial pipeline:              Sample_A/outs/spatial
- Run summary CSV:                          Sample_A/outs/metrics_summary.csv
- BAM:                                      Sample_A/outs/possorted_genome_bam.bam
- BAM index:                                Sample_A/outs/possorted_genome_bam.bam.bai
- Filtered feature-barcode matrices MEX:    Sample_A/outs/filtered_feature_bc_matrix
- Filtered feature-barcode matrices HDF5:   Sample_A/outs/filtered_feature_bc_matrix.h5
- Unfiltered feature-barcode matrices MEX:  Sample_A/outs/raw_feature_bc_matrix
- Unfiltered feature-barcode matrices HDF5: Sample_A/outs/raw_feature_bc_matrix.h5
- Secondary analysis output CSV:            Sample_A/outs/analysis
- Per-molecule read information:            Sample_A/outs/molecule_info.h5
- Loupe Browser file:                       Sample_A/outs/cloupe.cloupe
- Spatial Enrichment using Moran''s I file:  Sample_A/outs/spatial_enrichment.csv

Pipestance completed successfully!

输出的结果存放在一个以样本名为名称的文件夹内（如 Sample_A), 主要结果均位于 outs 文件夹内。

└── outs
    ├── analysis
    │   ├── clustering # 聚类文件夹，图聚类，k-means 聚类
    │   ├── diffexp # 差异分析
    │   ├── pca # pca 线性降维
    │   ├── tsne # tsne 非线性降维
    │   └── umap # umap 非线性降维
    ├── cloupe.cloupe # 用于 Loupe Browser 可视化和分析文件
    ├── filtered_feature_bc_matrix # 组织相关 spot 矩阵信息
    │   ├── barcodes.tsv.gz
    │   ├── features.tsv.gz
    │   └── matrix.mtx.gz
    ├── filtered_feature_bc_matrix.h5 # 过滤后的 barcode 信息 HDF5 格式
    ├── metrics_summary.csv # 结果汇总 csv 文件
    ├── molecule_info.h5  # UMI 信息，spaceranger aggr aggregate 数据的时候会用到的文件
    ├── possorted_genome_bam.bam
    ├── possorted_genome_bam.bam.bai
    ├── raw_feature_bc_matrix # 所有 spot 完整矩阵信息（组织和背景）
    │   ├── barcodes.tsv.gz
    │   ├── features.tsv.gz
    │   └── matrix.mtx.gz
    ├── raw_feature_bc_matrix.h5 # 所有 spot 信息 HDF5 格式
    ├── spatial # 空间信息相关信息:这些文件是用户提供的原始全分辨率 brightfield 图像的下采样版本。下采样是通过 box 滤波实现的，它对全分辨率图像中像素块的 RGB 值进行平均，得到下采样图像中一个像素点的 RGB 值。
    │   ├── aligned_fiducials.jpg # 这个图像的尺寸是 tissue_hires_image.png。由基准对齐算法发现的基准点用红色高亮显示。此文件对于验证基准对齐是否成功非常有用。
    │   ├── detected_tissue_image.jpg # 检测到的组织图像
    │   ├── scalefactors_json.json # 存放将各个图片文件与原始无损图片比例系数的文件
    │   ├── tissue_hires_image.png # 图像的最大尺寸为 2000 像素
    │   ├── tissue_lowres_image.png # 图像的最大尺寸为 600 像素
    │   └── tissue_positions_list.csv  # spot 坐标信息
    └── web_summary.html 结果汇总 html 文件

$ cat scalefactors_json.json
# {"spot_diameter_fullres": 71.45230555934738, "tissue_hires_scalef": 0.19729704, "fiducial_diameter_fullres": 115.4229551343304, "tissue_lowres_scalef": 0.05918911}

scalefactors_json.json：存放各个图片文件与全分辨率图片比例系数的文件。

1. 基准点直径在全分辨率图中所占的像素点目;

2. 全分辨率图与 high res 图片的缩放比例系数;

3. spot 直径在原片中所占的像素点目;

4. 全分辨率图与 low res 图片的缩放比例系数。

$ head -n 2 tissue_positions_list.csv 
# ACGCCTGACACGCGCT-1,1,0,0,1572,1178
# TACCGATCCAACACTT-1,1,1,1,1667,1234

tissue_positions_list.csv: spot坐标信息。每一行为一个spot的信息。

1. barcode: 与芯片相关的条码序列;

2. in_tissue: 二进制，表示这个点是在组织的内部(1)还是外部(0);

3. array_row: 数组中从 0 到 77 的点的行坐标。数组有 78 行;

4. array_col: 为了表示点的橙色板条箱排列，这个列索引对偶数行使用从 0 到 126 的偶数，对奇数行使用从 1 到 127 的奇数。注意，每一行(偶数或奇数)有 64 个位置;

5. pxl_col_in_fullres: 全分辨率图片中点中心的列像素坐标;

6. pxl_row_in_fullres: 全分辨率图片中点中心的行像素坐标。

web_summary.html：网页版报告对结果进行汇总并进行可视化：分为 Summary 和 Analysis 两个部分。

Summary 部分包含如下几个部分

1. 代表性指标和测序质量

对样本中的 spot 和表达的基因数评估，同时给出了测序获得的总的双端 reads 数、有效 barcodes 的百分比、有效 UMI 百分比、测序饱和度以及各种 Q30 百分比信息。如下图所示：

2. 比对质量

统计 reads 比对到参考基因组中的百分比，同时给出唯一比对到参考基因组中的 reads 百分比，唯一比对到基因间区、内含子区、外显子区和转录本区的 reads 百分比等。

3. Spots 信息

有关 spot 信息的统计，包括在组织下的 spot 的 reads 比例、每个 spot 的平均的 reads 数、在组织下的每个 spot 的平均 reads 数、在组织下的 spot 对应的 barcode 检测到的中位基因数、检测到的总的基因数和每个 spot 的中位 UMI counts 数。

4. 样本信息

有关样本分析的基本信息: 样本名、样本描述、文库类型、玻片序列号和区域、参考基因组信息和spaceranger 版本信息。

Analysis 部分展示如下结果

1. UMIs 分布展示

左图：图像上 UMI 的分布，具有越多的 UMI 的 spot 含有更多的 mRNA。

右图：t-SNE 降维可视化后的 UMI 的分布。在这个空间中，彼此靠近的成对 spots 比彼此远离的 spots 具有更相似的基因表达谱。

鼠标放置到图像上会现在对应的位置信息和对应 spot 上的 UMI count 数。从这个图我们可以判断 UMI 主要分布在组织的哪些区域、哪些区域没有捕获到 mRNA 或捕获的 mRNA 特别少。

2. 降维聚类结果展示

左图：cluster 在组织图像上的分布。

右图：t-SNE 降维可视化后的 cluster 的分布。

鼠标放置到图像上会现在对应的位置信息和对应 spot 上的 cluster 值和该 cluster 占总的 spot 的比例。

3. 基因差异表达分析

差异表达分析旨在为每个亚群找到在该亚群中相对于样本其余亚群表达更高的基因。这里对每个基因在每个亚群和样本其余亚群之间进行差异表达分析。表头的具体含义可点击右侧？进行了解。

4. 测序饱和度评估

对 reads 抽样，观察不同抽样条件下检测到的转录本数量占检测到的所有转录本的比例（测序饱和度)，并绘制如下曲线。测序饱和度是观察到的文库复杂度的一个度量，当所有转换的 mRNA 转录物都测序后接近 1.0（100%），即如果曲线末端区域平滑，说明测序接近饱和，再增加测序量，覆盖到的转录本数目也不会有明显增多。