分享是一种态度
单细胞RNA-seq分析介绍 单细胞RNA-seq的设计和方法 从原始数据到计数矩阵 差异分析前的准备工作 scRNA-seq——读入数据详解 scRNA-seq——质量控制 为什么需要Normalization和PCA分析 scRNA-seq聚类分析(一) scRNA-seq聚类分析(二) scRNA-seq Clustering (一) scRNA-seq Clustering (二) scRNA-seq Clustering quality control (一) scRNA-seq Clustering quality control (二)
现在,我们已经确定了所需的群集,可以继续进行标记识别,这将使我们能够验证某些群集的身份并帮助我们推测任何未知群集的身份。
我们的聚类分析产生了以下群集:
我们在聚类分析中有以下问题:
我们可以使用Seurat探索几种不同类型的标记识别,以获得这些问题的答案。每种都有自己的优点和缺点:
通常建议在评估单个样本组/条件时使用此类型的分析。通过 FindAllMarkers()
函数,我们将每个群集与所有其他群集进行比较,以识别潜在的标记基因。每个群集中的细胞被视为重复的,本质上是通过一些统计检验来执行差异表达分析。
注:默认值为Wilcoxon Rank Sum检验,但也有其他可用选项。
FindAllMarkers()
函数有三个重要参数,它们提供了确定基因是否为标记基因的阈值:
logfc.threshold
:相对于所有其他群集组合中的平均表达,群集中基因的平均表达的最小log2倍数变化。默认值为0.25。min.diff.pct
:群集中表达基因的细胞百分比与所有其他簇中表达基因的细胞百分比之和的最小百分比差异。min.pct
:只测试在两个群体中任何一个的细胞中检测到的最小部分的基因。旨在通过不测试那些很少表达的基因来加快功能。默认值为0.1。您可以根据需要的严格程度使用这些参数的任意组合。同样,默认情况下,此功能将返回给您同时显示正向和负向表达变化的基因。通常,我们添加一个参数only.pos
以选择仅保留积极的更改。查找每个集群标记的代码如下所示。我们不会运行此代码。
## DO NOT RUN THIS CODE ##
# Find markers for every cluster compared to all remaining cells, report only the positive ones
markers <- FindAllMarkers(object = seurat_integrated,
only.pos = TRUE,
logfc.threshold = 0.25)
注意:此命令可能要花很长时间才能运行,因为它正在针对所有其他细胞处理每个单独的群集。
因为我们的数据集中有代表不同条件的样本,所以我们最好的选择是找到保守的标记。此函数按样本组/条件在内部分离出细胞,然后针对所有其他群集(或第二个群集,如果指定,则为第二个群集)执行单个指定群集的差异基因表达测试。针对每种情况计算基因水平的p值,然后使用MetaDE R软件包中的meta分析方法进行跨组组合。
在开始标记鉴定之前,我们将明确设置默认测定,我们希望使用原始计数,而不是集成数据。
DefaultAssay(seurat_integrated) <- "RNA"
注意:虽然此函数的默认设置是从“RNA”插槽获取数据,但我们建议您运行上面的代码行,以绝对确保万一活动插槽在您的分析中的上游某处发生更改。原始计数和归一化计数存储在此槽中,用于查找标记的函数将自动提取原始计数。
函数 FindConservedMarkers()
具有以下结构:FindConservedMarkers()
语法:
FindConservedMarkers(seurat_integrated,
ident.1 = cluster,
grouping.var = "sample",
only.pos = TRUE,
min.diff.pct = 0.25,
min.pct = 0.25,
logfc.threshold = 0.25)
从中可以看到我们前面为 FindAllMarkers()
函数描述的一些参数;这是因为它在内部使用该函数首先在每个组中查找标记。在这里,我们列出了在使用 FindConservedMarkers()
时提供的一些附加参数:
ident.1
:该函数一次只评估一个集群;在这里您将指定感兴趣的集群。grouping.var:
元数据中的变量(列标题),它指定将细胞分隔成组对于我们的分析,我们将相当宽松,只使用大于0.25的log2倍数变化的阈值。我们还将指定仅返回每个群集的正标记。
在一个群集上进行测试,看看它是如何工作的:
cluster0_conserved_markers <- FindConservedMarkers(seurat_integrated,
ident.1 = 0,
grouping.var = "sample",
only.pos = TRUE,
logfc.threshold = 0.25)
FindConservedMarkers()
函数的输出是一个矩阵,其中包含按我们指定的群集的基因ID列出的假定标记的排序列表,以及相关的统计数据。请注意,为每个组(在我们的 Case,Ctrl和Stim)计算相同的统计信息集,最后两列对应于这两个组中的组合p值。
注意:由于每个细胞都被视为复制细胞,这将导致每个组内的p值膨胀!一个基因可能有一个极低的p值<1e-50,但这并不能证明它是一个高度可靠的标记基因。
在查看输出时,我们建议寻找在 pct.1
和 pct.2
之间表达差异较大且logFC较大的标记。例如,如果 pct.1
=0.90和 pct.2
=0.80,则可能没有那么令人兴奋的标记。然而,如果 pct.2
=0.1,那么更大的差异将会很有说服力。同样,我们感兴趣的是表达该标记的大多数细胞是否在我感兴趣的群集中。如果 pct.1
较低,例如0.3,则可能没有那么有意义。如上所述,这两个参数也是运行函数时可能包括的参数。如上所述,这两个参数也是运行函数时可能包括的参数。
添加带有基因注释信息的列可能会很有帮助。为此,将此文件(https://github.com/hbctraining/scRNA-seq/raw/master/data/annotation.csv)下载到您的数据文件夹。然后将其加载到R环境中:
annotations <- read.csv("data/annotation.csv")
注意:如果您有兴趣了解我们是如何获得此注释文件的,请查看链接(https://hbctraining.github.io/scRNA-seq/lessons/fetching_annotations.html)。
首先,我们将带有基因标识符的行名转换为自己的列。然后,我们将此注释文件与来自FindConservedMarkers()
的结果合并:
# Combine markers with gene descriptions
cluster0_ann_markers <- cluster0_conserved_markers %>%
rownames_to_column(var="gene") %>%
left_join(y = unique(annotations[, c("gene_name", "description")]),
by = c("gene" = "gene_name"))
View(cluster0_ann_markers)