转自Seurat_Satija 如何使用 Seurat 分析单细胞测序数据( Q&A
1. 关于 QC
Q:过滤线粒体基因是去除低质量 cell,核糖体基因用过滤吗?
A: Seurat 教程中没有提到用核糖体基因的 UMI 含量来过滤细胞,我看到的文献中基本上也不这么做。单细胞中表达量最高的基因一般是线粒体基因、核糖体基因以及 actin 和 hemoglobin基因。有人认为要去掉这些基因再做后续分析,也有人说不要去掉,没有一个肯定的答案。
2. 关于 NORMALIZATION
Q: 请问不同样本的测序深度用 NormalizeData 可以完全去除吗? 样本比较多的话,去除样本测序深度推荐用 SCT 还是标准的 Normalize 方法呢? / SCTransform 什么时间推荐使用?
A: 在用 NormalizeData 函数做 normalization 的时候, 每个细胞中所有基因的表达量除以该细胞中检测到的所有 count 数, 这样会使每个细胞所有基因的表达量之和都变得一样,也就是说细胞之间的测序深度都保持一致。但是这样做是否能够完全去除测序深度对基因表达变异的影响,还是要打个问号。 根据 Seurat 官网上的介绍: Even after standard log-normalization,variation in sequencing depth is still a confounding factor, and this effect can subtly influence higher PCs。 因此, Seurat 又提供了 SCTransform 的方法,并声称: In sctransform, this effect is substantially mitigated. This means that higher PCs are more likely to represent subtle, but biologically relevant, sources of heterogeneity – so including them may improve downstream analysis. 除了可以保留更多的生物学变异, SCTransform 对参数的依赖度更小,产生的结果更稳定。 根据这些优点, Seurat 推荐使用 SCTransform 做 normalization。
Q: 使用 SCTransform 的時候, VlnPlot 所呈現的數值似乎會偏小,視覺上似乎較不明顯,請問是否會影響分析結果?
A: 这里所说的数值偏小,是指同一个基因,用 SCTransform 比用标准的 normalization 方法得到基因表达量偏小吗? 我觉得如果整体上基因的数值都偏大或偏小,应该不会对后续分析有太大影响。不过,我还是建议你用两种 normalization 方法都跑一遍完整的流程,看看结果会不会有很大差异。
Q: Normalization 时,为什么要用 RPM 这种方式而不是像 DEseq2 计算 sizefactor 呢?
findmarker 时,如果用 negativebinomial 这个模式,是用 count 数据,它会自己做
normalization 吗?
A: 在用 FinderMarkers 函数做差异表达分析的时候,如果选择的是 DEseq2,函数内部会使用estimateSizeFactor 函数计算 sizefactor。 如果你使用的是 negbinom 检验, FindMarkers 内部调用的是 MASS::glm.nb 函数,但是我没有在 Seurat 源代码中看到在做这个检验之前有normalization 的步骤。
3. 关 于 特 征 选 择
Q: Findvariable 函数中 features 为什么要设置 2000?可以变吗(彩环 田) / 想请问您对挑选特征基因有什么看法? 预设 2000 HVGs 适用于大型资料吗?
A: FindVariableFeatures 函数有 3 种选择高表达变异基因的方法,可以通过 selection.method参数来选择,它们分别是: vst(默认值), mean.var.plot 和 dispersion。 nfeatures 参数的默认值是 2000,可以改变。如果 selection.method 参数选择的是 mean.var.plot,就不需要人为规定高表达变异基因的数目,算法会自动选择合适的数目。 我建议使用完 FindVariableFeatures 函数后,用 VariableFeaturePlot 对这些高表达变异基因再做个可视化,看看使用默认值 2000 会不会有问题。
Q:请问整合时,如果每个待整合样本中的基因个数小于 2000,参数 feature 该如何设置呢?
A: FindVariableFeatures 函数的 nfeatures 参数可以选择各个样本中基因数最小的那个值,但是在此之前,建议你先确认一下样本中的基因数比较少的原因,样本的质量是否过关, QC 过程中有没有发现问题, subset 函数筛选的时候有没有出错。
4. 关于细胞周期
Q: Cell Ranger 可以实现去除细胞周期影响吗?
A: Cell Ranger 不能去除细胞周期影响,建议你用 Seurat 读入 Cell Ranger 的输出文件,用CellCycleScoring 函数和 ScaleData 函数去除细胞周期影响。
Q: 请问能否再解释下 两种细胞周期去除方法的区别和意义?请问 scaledata 在去除细胞周期影响时候,您刚才讲过用 var.to.gress(S-G2M),或者 var.to.gress( S, G2M) ,这两种在意义上是什么区别?
A: 我们利用 ScaleData 函数的 vars.to.regress 参数来消除细胞周期带来的变异。这里可以用到两种方法,一种方法是同时 regress out S.score 和 G2M.score,另一种方法是 regress out 两个score 的差。从 PCA 图上可以看到,第一种方法完全去除了细胞周期的影响,而第二种方法则消除了 S 和 G2M 期的差别,为什么要这么做呢?因为对于某些样本来说,有的细胞类型处于休眠期,有的细胞类型处于增殖期,如果把所有细胞周期的影响都去除掉,会影响这两类细胞的鉴定。这时候,我们只要去除增值期内部细胞周期的差异,保留增值期和休眠期的差别就可以了。
Q: findmarker 使用的是标准化的 slot, data 么?如果是,那之前 scale 中去除的细胞周期影响,是否还有意义?
A: 如果用细胞周期 marker 基因做 PCA 分析,发现细胞可以明显按照细胞周期分为不同的group,说明细胞周期对这个数据中的基因表达量有比较大的影响,用 ScaleData 函数去regress out 细胞周期 score 是十分必要的。 Findermaker 函数使用的是 data 值,不是scale.data。如果在做差异表达分析的时候也考虑细胞周期的影响,可以使用 metadata 中的Phase 变量,先筛选处于相同细胞周期的细胞再做差异表达分析。
Q: scale 里对细胞周期的矫正,和 integrate 的整合矫正有哪些差异?
A: 做细胞周期矫正,首先要有细胞周期 marker 基因,然后用 CellCycleScoring 函数来计算细胞周期 score,然后用 ScaleData 函数去 regress out 细胞周期 score。 IntegrateData 主要是矫正多个样本间的批次效应,例如整合多个来自不同单细胞测序技术的数据
5. 关于降维分析
Q: ElbowPlot 的拐点是指这个点后面基本持平了么? / 請問如何判斷 elbow 拐點呢?為何不是在 3~4,而是在 9~10/ PC 拐点如何界定?
A: ElbowPlot 作图后,用肉眼观察,当曲线进入平台期,没有明显继续下降的趋势就可以了。ElbowPlot 默认显示 20 个 PC,如果发现 20 个 PC 还没有到达平台期,可以修改 ElbowPlot 函数的 ndims 参数。 RunPCA 函数默认计算 50 个 PC,如果想看更多的 PC,可以修改 RunPCA 的npcs 参数。
Q:请问 PCA 那页,是细胞和基因都做过 PCA 么?
A: 在 Seurat 中,我们是对细胞做 PCA 分析,也就是对细胞做降维分析,比如从 2000 维(基因)降到 50 维( PC)。一般情况下,我们不会对基因做 PCA 分析。
Q: PCA 中的 embedding 可以再解释下么?
A:比如,原本每个细胞有 2000 个基因,每个基因在各个细胞中都有表达量,做了 PCA 分析后,每个细胞有 50 个 PC,每个 PC 在各个细胞中都有 embedding 值。为了便于理解,你可以把 PC 看成是整合了多个基因的 meta-gene,把 cell embedding 看成是这些 meta-gene 的表达量。
Q: 可不可以讲下 umap 和 tSNE 的主要原理?
A: UMAP 和 tSNE 都是非线性的降维方法,是把在高维空间中的点(我们人类肉眼无法观察)投射到二维平面上。这些高维空间中的点是有一定结构的,有些点彼此之间距离更近,有些点距离更远。我们在降维的过程要尽量保留这种结构,让在高维空间中更近的点,在二维平面上也更近,如果降维分析打乱了数据原有的结构,就毫无意义了。 tSNE 的目标是在降维的过程中尽量保留数据的局部结构,而 UMAP 除了能保留数据的局部结构,还能更好的保留全局结构。如果想更深入的了解它们的原理,需要一定的数学基础,我推荐两篇解释的相对简单一些的文章供你参考:
https://opentsne.readthedocs.io/en/latest/tsne_algorithm.html,
How UMAP Works - umap 0.4 documentation
Q: 我想问下数据降维为什么先做 PCA 然后再 UMAP 的降维 这俩不都是降维吗? 不能直接
UMAP?
A: PCA 分析和 UMAP 分析的目的不同。 PCA 分析把高表达变异基因整合成若干个 PC。这样做一方面降维,一方面降低噪音。做聚类分析和 UMAP 分析的时候,使用的不是基因表达量,而是由 PCA 分析得到的 PC 的值。 UMAP 的主要作用是把高维度中距离接近的点映射到二维平面上距离相近的位置上。 UMAP 返回的结果只有两个维度( UMAP_1 和 UMAP_2),只能用来做数据可视化。
6. 关 于 聚 类 分 析
Q: cluster ID 为何从 0 开始?是否可以从 1 开始?从 0 开始是因为有什么‘算法’方面特殊的优势么?
A: 默认 clusterID 是从 0 开始,没有特殊的优势,可以使用 RenameIdents 函数重命名clusterID。
Q: 鉴定出来的细胞类群用 UMAP 绘图后发现聚类效果很差,大量 cluster 相互交叉, 同一类群细胞弥散分布是什么原因?聚类可视化和类群划分不匹配。
A: 首先要确保聚类和 UMAP 用的 dims 参数是一致的。可以适当降低一下 FindClusters 函数的resolution 参数,减少 cluster 数目,看看能不能把相互交叉的 cluster 聚成一个 cluster。 还可以尝试 FindClusters 函数中不同的 algorithm 参数,看看聚类效果会不会改进。
Q: 请问分群过程中,是否有设置邻居个数的函数和参数?
A: FinderNeighbors 函数的 k.param 参数是用来设置 KNN 算法中最近邻的个数。
Q: 请问 finder cluster 只能使用 SNN 进行聚类么?可以有其他选择吗? / seurat 的聚类方式除了 KNN 外还有其他的选择吗?
A: Seurat 的聚类方法是基于 SNN 图和 Louvain 或 SLM 算法, FindNeighbors 函数返回的SNN 图是在 KNN 图的基础上得来的,不支持其他方法。
7. 关 于 数 据 整 合
Q:整合多样本分析时,有些样本细胞数量过少,可能不足 200,这样在整合数据时默认参数就会报错,有什么好的建议吗?亚群细分时也经常出现这个问题,及个别样本细胞数量过少,不能进行批次校正/ 在对整合样本中的某一细胞亚群进行再次分群时,如果有些样本中这一亚群的细胞比较少,怎么进行整合分析?
A:在不知道具体报错信息的情况下,我无法判断问题的原因是由于细胞数目少还是来自样本本身。 我的建议是: 首先要确定不足 200 个细胞的样本是否都会遇到这个问题。如果你在其他样本中随机抽取不到 200 个细胞,去做整合,也会报同样的错误吗? 如果给这些不足 200 个细胞的样本逐渐增加细胞数,达到多少细胞数的时候,就不会报错?
Q:请问整合时,报这样的错该如何解决? ‘Error in irlba(A = t(x = object), nv =npcs, ...) :max(nu, nv) must be positive
A:我看不到你的代码和数据,所以不好判断是什么问题,有可能是输入数据有问题,你可以检查一下输入的数据是不是都做了相同的预处理,包括 normalization 的方法和高变基因的数目等等,再看看 FindIntegrationAnchors 和 IntegrateData 函数中选择的 dims 参数是否一致。如果你整合的样本比较多,可以先拿两个样本来整合,如果成功,再继续增加样本,直到报错,就可以定位是哪个数据出了问题。
Q:在整合多个样品的单细胞数据时系统报错,说超过内存怎么办?我们的样品大约有 20 多个/ 在整合多个样本数据时,由于样本量过大,导致超过稀疏矩阵上限,无法进行 integrateData步骤,应如何解决?
A:我也会经常遇到数据量超过内存上限的问题,解决办法就是换一个有更大内存的服务器,此外,在做数据分析的过程中,经常查看当前环境中的数据占用了多大内存,及时清理一些中间变量,删除之前最好将一些重要的变量保存在 RDS 文件里。如果这些方法都不行,就只好对样本中的细胞随机抽样,选取一部分细胞来做整合。
Q:样本之间如果没有 batch effect,可以用直接 merge 合并吗? / 我们可不可以直接把多个sample 直接合并不用 intergation 函数如果没看出有 batch effect/ 关于样本整合, merge 函数可以用吗?
A:关于样本整合,建议使用 Seurat 提供的 FindIntegrationAnchors 和 IntegrateData 函数。merge 函数只是合并多个 Seurat 对象,无法消除批次效应。如果你的数据没有批次效应,merge 也是可以用的。
Q:多个样本整合后的细胞分群结果是不是可能和单个样本做细胞分群结果会不一致?
A: 这是很有可能的。比如一个样本中某种细胞的含量很少,它们有可能和其他类型的细胞聚在一起,但是另一个样本中该类细胞的含量比较多,当两个样本整合在一起后, 两个样本中的这类细胞应该更倾向于聚在一起,这样就和单个样本做细胞分群的结果不一致。这只是我的假设,我没有用实际数据去验证过这种不一致的程度有多大,和哪些因素有关。你如果有这样的数据可以尝试分析一下。
Q: 对于不同时间,同样用 10X 处理的单细胞转录组样本,但是测序饱和度不同,必须用integration 方法进行批次效应处理么?还是说直接合并就可以了? / 请问多个发育时期的样本需要整合流程吗?
A: 首先要判断样本之间是否有批次效应。一般是做一个 tSNE 或 UMAP 图,如果来自同一个样本的细胞都聚在一起,来自不同的样本的细胞都分的很开,说明批次效应明显,需要用integration 的方法进行数据整合。整合完成后再做一个 tSNE 或 UMAP 图,查看批次效应是否已经去除。如果没有批次效应,直接 merge 也是可以的。 如果在 tSNE 或 UMAP 图上不容易判断是否有批次效应,建议两种方法都试一下,看看最后的结果有没有很大差异。
Q: 多样本整合是按样本来整合吗? 还是按不同样本类型来整合?
A: Seurat 可以整合不同的样本类型,比如 scRNA-seq 和 CITE-seq 的整合, 以及 scRNA-seq和 scATAC-seq 的整合。可以参考 2019 年的那篇 Cell 文章。
Q: CITE-seq 用 seurat 分析如何将多个样本合并呢?
A: 将 CITE-seq 的 ADT 数据存放在 Assay 对象中,并设置为默认 assay,然后用和 scRNA-seq一样的方法整合数据。
Q: 整合会不会把样本间的部分差异当成批次效应给去除掉?
A:我觉得是有可能的,这种“矫枉过正”的程度有多大是和样本之间本身的差异有关的。可以设想一下,如果把 PBMC 的样本和神经细胞的样本放在一起整合会是什么效果。这两个样本的生物学差异应该大于批次效应,这时候去批次效应就有可能把样本间的部分差异也去除掉了。 这只是我的推测,并没有用数据验证过,仅供参考。
Q: 如果不进行 intergration ,可以直接比较 control 和疾病组的差异吗?整合的优势是统一标准化吗?
A:如果是 bulk RNA-seq,可以直接比较 control 和疾病组的差异。 但是 bulk RNA-seq 中基因的表达量是基因在所有细胞中表达量的平均值,没有考虑到样本中细胞的异质性,某个基因在某种细胞类型中的差异表达有可能被样本中细胞类型的差异所掩盖。因此, 我们要利用单细胞RNA-seq,找到两个样本中相同类型的细胞后,再寻找这类细胞在 control 和疾病组的差异表达基因。具体来说,我们利用 Seurat 整合完 control 和疾病组样本,并做聚类分析后,分别在control 和疾病组中对某个 cluster 做差异表达分析,寻找这个 cluster 的 marker 基因,然后,把在 control 和疾病组中都被鉴定为 marker 基因的基因, 作为该 cluster 保守的 marker 基因,并利用这些基因对细胞类型进行注释。最后,针对同一类细胞,寻找在 control 和疾病组中发生变化的基因。
Q: 如果我想对某种细胞分亚群,是用经 PCA 和的 seurat 对象取出这部分细胞聚类,还是需要重新整合,或者重新 PCA?
A: 在分亚群的时候,很多人都会遇到类似的问题。在 GitHub 上也有很多讨论:
re_clustering in seurat v3 · Issue #1528 · satijalab/seurat
How to analysis the subset cells of already intergrated object? · Issue #1547 · satijalab/seurat
但是,目前还没有看到 Seurat 官方给出的最佳实践。我的建议是在分亚群之前,先在 UMAP 图上看看有没有 batch effect,如果没有,就不需要重新整合。
8. 关 于 差 异 表 达 分 析
Q: 老师您好,请问差异分析为什么使用 RNA assay ?而不是用 SCT assay 或者 Integratedassay? RNA assay 不是未经过标准化和批次校正的原始数据么? / 整合数据后进行 findmarker分析时是采用 integrated assay 还是 RNA assay?
A: 这个问题也困扰了我很久,后来我在 GitHub 上找到了官方的解释: As a general rule, wealways recommend performing DE on originally measured values - instead of on batchcorrected, imputed, etc. values. This ensures that the measurements that enter the DE test are indeed independent from each other, which is a requirement of any statistical DE test. Find DEGs after doing integration · Issue #1256 · satijalab/seurat
Q: 得到的差异基因和保守基因,用哪个作为 marker 基因好? / FindallMarkers 和findconservedmarkers 结果上有什么区别吗?
A: 如果只有一个样本,就用 cluster 中的差异表达基因做 marker 基因。如果做了数据整合,比如处理组和对照组,就要先用 FindConservedMarkers 函数找到保守的差异表达基因,用它们来作为 marker 基因对 cluster 做注释,注释完以后,再用 FindMarkers 函数比较这个 cluster 中的处理组细胞和对照组细胞之间的差异表达基因。
Q: 如果想从所有基因中筛选差异基因,有方法能实现么?
A: 在用 FindMarkers 函数找差异表达基因的时候, 可以通过 logfc.threshold(默认值 0.25),min.pct(默认值 0.1)和 min.diff.pct(默认值-Inf)三个参数来预筛选基因,这样可以加快计算速度,如果想从所有基因中筛选差异基因,就把 logfc.threshold 设为 0, min.pct 设为 0。
9. 关于细胞注释
Q: cluster 名称的命名主要参考哪几方面的信息(如 cluster 的 top makers)?
A:对 cluster 进行注释,主要是看 FindMarkers 或 FindAllMarkers 函数找到的 marker 基因。利用已有的知识或者数据库,将 cluster 命名为特定的细胞类型。有几个细胞 marker 基因的数据库供参考:
CellMarker: CellMarker
PanglaoDB: A Single Cell Sequencing Resource For Gene Expression Data
Q: 您可以再详细讲一下如何自动的识别细胞类型吗?
A:首先是用 FindIntegrationAnchors 函数找到数据集间的锚点,然后用 IntegrateData 函数整合多个数据集。这个数据集可以用来研究处理组和对照组的变化。如果这个数据集有细胞类型的注释信息,就可以拿来作为参考数据集,对查询数据集进行细胞类型的注释。方法也很简单,先使用 FindTransferAnchors 函数,找到参考数据集和查询数据集的锚点,然后使用TransferData 函数和参考数据集的标签对查询数据集中的细胞进行注释,每个细胞都有一个预测的细胞类型以及打分。这个预测信息可以用 AddMetaData 函数添加到查询数据集的metadata 中。
Q:请问癌症单细胞分群后注释,怎么区分恶性细胞和正常细胞群?
A: 我所知道的方法有两种,一种是看 marker 基因,另一种是用基因表达数据推断 CNV。 给你两篇文献供参考:
(1) Single-Cell Transcriptomic Analysis of Primary and Metastatic Tumor Ecosystems in Head and Neck Cancer.
https://www.cell.com/cell/fulltext/S0092-8674(17)31270-9
(2) Single-cell RNA sequencing demonstrates the molecular and cellular reprogramming ofmetastatic lung adenocarcinoma.
http://www.nature.com/articles/s41
原文:如何使用 Seurat 分析单细胞测序数据( Q&A)-上 - 知乎 (zhihu.com)