单细胞分析:质控实操(五)

  • Post author:
  • Post category:其他




1. 学习目标




  1. 构建质量控制指标并评估数据质量


  2. 适当的应用过滤器去除低质量的细胞



2. 过滤目标



  1. 过滤数据以仅包含高质量的真实细胞,以便在对细胞进行聚类时更容易识别不同的细胞类型
  2. 对一些不合格样品的数据进行检查,试图查询其不合格的原因



3. 挑战



  1. 从少量复杂的细胞中描绘出质量较差的细胞
  2. 选择合适的过滤阈值,以便在不去除生物学相关细胞类型的情况下保留高质量的细胞



4. 质量标准



当数据加载到

Seurat

并创建初始对象时,会为计数矩阵中的每个单元组装一些基本元数据。要仔细查看此元数据,查看存储在

merge_seurat

对象的

meta.data

中的数据帧:

# 查看元数据
View(merged_seurat@meta.data)  # 具体介绍见质控准备章节
meta.data
meta.data

为了可视化质量控制分析情况,需要计算一些额外的指标。这些包括:


  • number of genes detected per UMI:

    这个指标让了解数据集的复杂性(每个

    UMI

    检测到越多基因,数据越复杂)

  • mitochondrial ratio:

    该指标将提供来自线粒体基因的细胞读数百分比



5. Novelty score



这个值很容易计算,取每个细胞检测到的基因数量的

log10

和每个细胞的

UMI

数量的

log10

,然后将

log10

的基因数量除以

UMI



log10

数量。

# 将每个单元格的每个 UMI 的基因数添加到元数据
merged_seurat$log10GenesPerUMI <- log10(merged_seurat$nFeature_RNA) / log10(merged_seurat$nCount_RNA)



6. 线粒体率




Seurat

有一个方便的功能,可以计算映射到线粒体基因的转录本比例。

PercentageFeatureSet()

函数接受一个模式参数,并在数据集中的所有基因标识符中搜索该模式。由于正在寻找线粒体基因,因此搜以“MT-”模式开头的任何基因标识符。对于每个细胞,该函数获取属于“Mt-”集的所有基因(特征)的计数总和,然后除以所有基因(特征)的计数总和。该值乘以 100 以获得百分比值。

# 计算百分比
merged_seurat$mitoRatio <- PercentageFeatureSet(object = merged_seurat, pattern = "^MT-")  # 该模式 ^MT- 应该根据自己的数据集进行修改
merged_seurat$mitoRatio <- merged_seurat@meta.data$mitoRatio / 100

现在已经具备了评估数据所需的质量指标。但是,希望在元数据中包含一些有用的附加信息,包括单元 ID 和条件信息。首先通过从

Seurat

对象中提取

meta.data

来创建元数据:

# 创建元数据
metadata <- merged_seurat@meta.data

# 将 Cell ID 添加到元数据
metadata$cells <- rownames(metadata)

# 创建样本列
metadata$sample <- NA
metadata$sample[which(str_detect(metadata$cells, "^ctrl_"))] <- "ctrl"
metadata$sample[which(str_detect(metadata$cells, "^stim_"))] <- "stim"

# 重命名列
metadata <- metadata %>% dplyr::rename(seq_folder = orig.ident,
                                       nUMI = nCount_RNA, 
                                       nGene = nFeature_RNA)
# 最终结果如下图
metadata
metadata
# 将更新的元数据保存到` Seurat` 对象
merged_seurat@meta.data <- metadata

# 保存为`.RData`
save(merged_seurat, file="data/merged_filtered_seurat.RData")
# 结果如下
文件结构
文件结构



7. 质量评估指标



下面将评估以下各种指标,然后决定哪些

cells

质量低,应从分析中删除


  • Cell counts

细胞计数由检测到的独特细胞条形码的数量决定。对于本实验,

预计在 12,000 -13,000 个细胞之间

在理想情况下,会期望唯一细胞条形码的数量与加载的细胞数量相对应。然而,情况并非如此,因为细胞的捕获率只是加载的一部分。例如,与 50-60% 之间的

10X

相比,

inDrops

细胞捕获效率更高(70-80%)。

细胞数量也可能因

protocol

而异,产生的细胞数量远高于加载的数量。例如,在

inDrops protocol

期间,细胞条形码存在于水凝胶中,这些水凝胶与单个细胞和裂解/反应混合物一起封装在液滴中。虽然每个水凝胶都应该有一个与之相关的细胞条形码,但有时水凝胶可以有多个细胞条形码。同样,使用

10X protocol

时,有可能仅在乳液液滴 (GEM) 中获得带条形码的珠子,而没有实际的细胞。除了死亡细胞的存在之外,这两者都可能导致比细胞更多的细胞条形码。

# 可视化每个样本的细胞计数
metadata %>% 
   ggplot(aes(x=sample, fill=sample)) + 
   geom_bar() +
   theme_classic() +
   theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1)) +
   theme(plot.title = element_text(hjust=0.5, face="bold")) +
   ggtitle("NCells")
细胞计数
细胞计数

每个样本超过 15,000 个细胞,这比预期的 12-13,000 个要多。很明显,可能存在一些垃圾“细胞”。


  • UMI counts per cell

每个

spot



UMI

计数通常应高于 500,这是预期的。如果**

UMI

计数在 500-1000 计数之间,是可用的**,但可能应该对细胞进行更深入的测序。

# 可视化  UMIs/transcripts per cell 数量
metadata %>% 
   ggplot(aes(color=sample, x=nUMI, fill= sample)) + 
   geom_density(alpha = 0.2) + 
   scale_x_log10() + 
   theme_classic() +
   ylab("Cell density") +
   geom_vline(xintercept = 500)
nUMI
nUMI

可以看到,两个样本中的大多数细胞都有 1000 或更高的

UMI


  • Genes detected per cell

对基因检测的期望与对

UMI

检测相似,尽管它可能比

UMI

低一点。对于高质量数据,比例直方图应包含一个代表被封装细胞的大峰。如果看到主峰左侧有一个小肩,或者细胞的双峰分布,这可能表明有一些问题。可能有一组单元由于某种原因失败了。也可能是存在生物学上不同类型的细胞。

# 通过直方图可视化每个细胞检测到的基因分布
metadata %>% 
   ggplot(aes(color=sample, x=nGene, fill= sample)) + 
   geom_density(alpha = 0.2) + 
   theme_classic() +
   scale_x_log10() + 
   geom_vline(xintercept = 300)
nGene
nGene

  • Complexity

    (novelty score)

可以使用一种名为

novelty score

的来评估每个细胞

RNA

种类的复杂程度。

novelty score

是通过

nGenes



nUMI

的比率来计算的。如果有很多捕获的转录本(高

nUMI

)并且在一个细胞中检测到少量基因,这可能意味着只捕获了少量基因,并且只是一遍又一遍地从这些较少数量的基因中测序转录本。这些低复杂性(低

novelty score

)的细胞可能代表特定的细胞类型(即缺乏典型转录组的红细胞),或者可能是由于人为因素或污染造成的。一般来说,

预计优质细胞的

novelty score

高于 0.80

# 通过可视化每个 UMI 检测到的基因来可视化基因表达的整体复杂性(novelty score分)
metadata %>%
   ggplot(aes(x=log10GenesPerUMI, color = sample, fill=sample)) +
   geom_density(alpha = 0.2) +
   theme_classic() +
   geom_vline(xintercept = 0.8)
image-20221012214231411
image-20221012214231411

  • Mitochondrial counts ratio

该指标可以确定是否存在大量来自死亡或垂死细胞的线粒体污染。将线粒体计数的劣质样本定义为超过 0.2 线粒体比率标记的细胞。

# 可视化每个细胞检测到的线粒体基因表达分布
metadata %>% 
   ggplot(aes(color=sample, x=mitoRatio, fill=sample)) + 
   geom_density(alpha = 0.2) + 
   scale_x_log10() + 
   theme_classic() +
   geom_vline(xintercept = 0.2)
image-20221012214308103
image-20221012214308103

  • Joint filtering effects

孤立地考虑这些

QC

指标中的任何一个都可能导致对信号的误解。例如,线粒体计数比例较高的细胞可能参与呼吸过程,并且可能是想要保留的细胞。同样,其他指标可以有其他生物学解释。执行

QC

时的一般经验法则是将单个指标的阈值设置为尽可能宽松,并始终考虑这些指标的联合影响。通过这种方式,可以降低过滤掉任何活细胞群的风险。

经常一起评估的两个指标是**

UMI

的数量和每个细胞检测到的基因数量**。在这里,绘制了基因数量与线粒体读数分数着色的

UMI

数量的关系。联合可视化计数和基因阈值并额外覆盖线粒体分数,得出每个细胞质量的总结图。

# 可视化检测到的基因与 UMI 数量之间的相关性,并确定是否存在大量基因/UMI 数量少的细胞
metadata %>% 
   ggplot(aes(x=nUMI, y=nGene, color=mitoRatio)) + 
   geom_point() + 
 scale_colour_gradient(low = "gray90", high = "black") +
   stat_smooth(method=lm) +
   scale_x_log10() + 
   scale_y_log10() + 
   theme_classic() +
   geom_vline(xintercept = 500) +
   geom_hline(yintercept = 250) +
   facet_wrap(~sample)
image-20221012214520416
image-20221012214520416

好的细胞通常会表现出每个细胞更多的基因和更多的

UMI

(图的右上象限)。质量差的细胞可能每个细胞的基因和

UMI

较低,并且对应于图左下象限中的数据点。通过该图,评估了线的斜率,以及该图右下象限中数据点的任何散布。这些细胞具有大量的

UMI

,但只有少数基因。这些可能是垂死的细胞,但也可能代表低复杂性细胞类型(即红细胞)的群体。

线粒体分数仅在很少(颜色较深的数据点)的特别低计数的细胞中较高。这可能表明其细胞质

mRNA

已通过破裂的膜泄漏出来的受损/垂死细胞,因此,只有位于线粒体中的 mRNA 仍然是保守的。可以从图中看到,这些细胞被计数和基因数阈值过滤掉了。



8. 过滤




  • Cell-level

    过滤

现在已经可视化了各种指标,可以决定要使用的阈值,这将导致删除低质量的单元格。前面提到的建议通常是一个粗略的指导,具体的实验需要告知选择的确切阈值。下面将使用以下阈值:


  • nUMI > 500


  • nGene > 250


  • log10GenesPerUMI > 0.8


  • mitoRatio < 0.2

为了过滤,将回到

Seurat

对象并使用

subset()

函数:

# 使用选定的阈值过滤掉低质量的细胞
filtered_seurat <- subset(x = merged_seurat, 
                         subset= (nUMI >= 500) & 
                           (nGene >= 250) & 
                           (log10GenesPerUMI > 0.80) & 
                           (mitoRatio < 0.20))

  • Gene-level

    过滤

在数据中,将有许多计数为零的基因。这些基因可以显着降低细胞的平均表达,因此将从数据中删除它们。首先确定每个细胞中哪些基因的计数为零:

# 提取计数
counts <- GetAssayData(object = filtered_seurat, slot = "counts")

# 输出一个逻辑矩阵,为每个基因指定每个细胞的计数是否超过零
nonzero <- counts > 0

现在,将按

novelty score

进行一些过滤。如果一个基因只在少数几个细胞中表达,那么它并不是特别有意义,因为它仍然会降低所有其他不表达它的细胞的平均值。选择只保留在 10 个或更多细胞中表达的基因细胞。通过使用此过滤器,将有效去除所有细胞中计数为零的基因。

# 对所有 TRUE 值求和,如果每个基因超过 10 个 TRUE 值,则返回 TRUE
keep_genes <- Matrix::rowSums(nonzero) >= 10

# 只保留那些在超过 10 个细胞中表达的基因
filtered_counts <- counts[keep_genes, ]

最后,获取这些过滤计数并创建一个新的 Seurat 对象以进行下游分析。

# 重新分配给过滤后的 Seurat 对象
filtered_seurat <- CreateSeuratObject(filtered_counts, meta.data = filtered_seurat@meta.data)



9. 重新评估



执行过滤后,建议回顾指控指标以确保数据符合预期并且有利于下游分析。

本文由

mdnice

多平台发布



版权声明:本文为swindler_ice原创文章,遵循 CC 4.0 BY-SA 版权协议,转载请附上原文出处链接和本声明。