单细胞数据处理小细节汇总

kuaidi.ping-jia.net  作者:佚名   更新日期:2024-07-09

1. Seurat对象查看当前的Assay

在进行了SCTransform操作后,矩阵默认会变成SCT矩阵,如果不加设置,后续的PCA等操作都是基于SCT矩阵。
修改DefaultAssay:

2. Seurat使用FindVariables找到高变基因后增删高变基因

3. 不同运行步骤中的文件来源和储存位置⚠️

⚠️:PCA的值是可以被覆盖的,使用三步法对矩阵进行标准化后进行PCA后再使用SCT矩阵进行标准化,PCA的矩阵变成了SCT的PCA矩阵,原有的PCA矩阵不会保留。后续的TSNE和UMAP降维图也和三步法不一样。

4. @data标准化矩阵 和 @scale.data 归一化矩阵 的区别
单细胞RNA 测序数据中,文库之间测序覆盖率的系统差异通常是由细胞间的cDNA 捕获或PCR 扩增效率方面的技术差异引起的,这归因于用最少的起始材料难以实现一致的文库制备。标准化旨在消除这些差异(例如长度,GC 含量),以使它们不干扰细胞之间表达谱的比较。这样可以确保在细胞群体中观察到的任何异质性或差异表达都是由生物学而不是技术偏倚引起的(批次矫正仅在批次之间发生,并且必须同时考虑技术偏差和生物学差异,标准化只需考虑技术差异)。
软件Seurat 提供了三种标准化的方法,分别为LogNormalize、CLR、RC,通常情况下我们采用LogNormalize 的方式进行标准化,计算公式为:log1p((Feature counts/total counts) ∗ scale. factor)

归一化的目的则是使特征具有相同的度量尺度

参考: Seurat的normalization和scaling

5. 关于有些细胞属于同一个cluster但是在umap或者tsne图上相聚较远的问题:UMAP和TSNE是各自的算法在PCA降维的基础上再进行非线形降维,在二维图上把其各自算法认为相近的细胞聚在一起。但是FindClusters输入的不是UMAP或TSNE降维的数据,而是FindNeighbors的数据,而FindNeighbors输入的数据是PCA降维数据,是用另外一种算法计算的细胞之间的距离。因此会出现有些细胞被认为是同一个cluster,但是在umap或者tsne图上相聚较远(尤其是一些散在的,脱离主群的细胞)

6. marker基因鉴定,查看marker基因的表达是使用RNA矩阵还是sct矩阵?
这是一个争议性问题,两个都可以,目前建议最好使用RNA矩阵。
sct的到的count并不是真实的基因表达值,而是通过scaledata倒推出来的,它是一个回归,运算之后的残差。

7. 关于FindAllMarks找到的基因
如下图,先看cluster0的Marker基因:cluster0的差异基因是cluster0的细胞和剩下的所有的cluster合在一起的细胞做对比得到的。pct.1是这个基因在cluster0中的表达比例(S100A8在cluster0的细胞中的表达比例是100%),pct.2是这个基因在除了cluster0以外的所有细胞中的表达比例(S100A8在除cluster0以外的细胞中的表达比例是51.2%)。avg_log2FC是表达差异倍数,p_val_adj是校正后的p值。

8. 在提取Marker基因时比较好的办法:因为单细胞矩阵算出来的结果,p_val_adj==0的有很多,所以可以先把p_val_adj==0先提出来。再把p_val_adj<0.01的按差异倍数靠前的20/30/500...(按需要)个基因提出来,然后把两个矩阵合在一起(取交集)用来做细胞鉴定。(结合p值和fc来做筛选效果更好)

⚠️:提取没有核糖体和线粒体的marker基因更好。(这些基因对鉴定没有帮助)

有些基因比如Foxp3,对细胞鉴定很重要,但常常在筛选Marker基因的时候筛选不出来。 非负矩阵分解 可能更好。
参考: 过滤线粒体核糖体基因

9. 提取亚群

⚠️ 新提取的亚群需要重新进行降维聚类 (和大群相比,标记基因发生了变化),并重新寻找marker基因,重新分群,注释。❗️subset提取子集后,不同样本间批次校正的信息也被去除,需要重新进行批次矫正

参考: Seurat取子集时会用到的函数和方法 ⚠️⚠️⚠️

10. 取子集后如何去除空子集(还存在这个level,但里面包含的细胞为0,如何去除) as.factor(as.character())

11. 双细胞的预测和去除如DoubletFinder建议单样本进行,不建议双样本一起预测。除此之外,其他步骤都可以多样本一起做,质控的时候也可以多样本一起做,但是建议每个样本都单独看一看。

12. 单细胞多样本整合:merge();多样本拆分:SplitObject()

13. 在做多组数据整合,每个组又有多个样本的时候,最好把单独的每个样本当成批次,而不是把不同的组当成批次。

14. 多核运算

参考: 单细胞数据分析中future包的使用

15. pbmc3k.final@commands$FindClusters 可以查看FindClusters运行时间和记录。Seurat是记录其分析过程的,也可查看command下其他操作

16. 关于质控标准:同一组织的最好用同样的标准,不同组织的可以不一样。不同组织线粒体含量等可能存在差异。

17. 可视化的方法总结
参考: https://www.jianshu.com/p/0d1e2c7d21a4

18. circos图绘制

19. 单细胞数据思维导图,有利于查看单细胞数据格式。
https://www.jianshu.com/p/7560f4fd0d77

20. 对于旧版本Seurat对象的更新
scRNA <- UpdateSeuratObject(scRNA)
UpdateSeuratObject {SeuratObject} :Update old Seurat object to accommodate new features

21. 对有些操作需要用到python设置的情况

22. 单细胞数据做pooling的好处:可以尽量的降低dropout的问题。(dropout就是矩阵中的zero,这些zero实际上并不是0,而是每个液滴里面起始反应量太低了。而一般的反转录效率只能到30%左右,70%的转录本实际上在反转录那一步是被丢掉的,这是单细胞测序一个比较大的问题)。
但是一旦做了pooling,你必须要证明pooling对结果是没有影响的(下图的右面三个图)。

23. Seurat的VlnPlot中的combine参数,在如下画三个基因的情况下,设置成T就画一张图,设置成False,会将三个基因各画一张图。

24. rev()这一步是将横坐标的基因反过来排序

这两个画出来的图横坐标基因的顺序是相反的(见NicheNet)

25. 堆叠小提琴图的绘制
完成这个需求有以下几种实现方法:
1. Seurat包直接就可以实现(stack = T)
2. 通过Seuart->scanpy来实现,第一张是Seurat包VlnPlot函数画的图,第二张是scanpy中stacked_violin函数画的图,那么现在问题就变成为Seurat对象到scanpy对象的转换
3. 用R原生函数实现StackedVlnPlot的方法
4. 使用基于scanpy包衍生的scanyuan包来说实现
5. 使用R包MySeuratWrappers来实现

最简单的方法1如下:

如果不设置level,会按字母顺序排列,case会自动排在con前面。

使用Seurat的 RenameIdents 函数也可以

x: table
margin: a vector giving the margins to split by. E.g., for a matrix 1 indicates rows, 2 indicates columns, c(1, 2) indicates rows and columns. When x has named dimnames, it can be a character vector selecting dimension names.

得到的HC_1样本的orig.ident默认是样本名中第一个_号的前一部分。所以要保证矩阵的列名是 样本名_细胞barcode 这样的格式。
如果有多个分组,例如两个样本矩阵中细胞分别命名为HC_1_barcode,HC_2_barcode,在直接通过如下方法得到两个Seurat对象,再对其进行merge之后,两个样本会被合并成一个。也就是样本信息只保留了第一个_号之前的HC,没有保留_号之后的1和2。

为了避免这种情况,可以在构建Seurat对象时通过参数进行设置

⚠️PC数的选择:Seurat官网提供的三种方法只能给出PC数的粗略范围,选择不同PC数目,细胞聚类效果差别较大,因此,需要一个更具体的PC数目。作者提出一个确定PC阈值的三个标准:

一般先选默认分辨率(0.8),大概可能会分出十几个群。因为最终都是要注释到每一个barcode,所以首先可以看大类marker的分布(不受分辨率影响),可以根据marker基因的分布来调整分辨率。是否需要精细的分群得看精细的分群对研究有没有决定作用,还有很重要的一点是 看分出的各个cluster在Findallmarkers给出的结果中marker的热图是不是能明显分开 。精细划分的细胞本来就很类似,如果有部分小群的热图明显分不开或者非常类似,就可以考虑把分辨率调小。

这实际上是没有必要的必须保持一致的。下游的都是用pca之后的,pca是为了压缩数据。
umap和tsne是为了可视化(仅仅是可视化),但是FindNeighbor是计算细胞间距离矩阵。找类群数目和可视化可以说没有关系。

map函数:
R语言循环第三境界:purrr包map函数!
浅析R语言中map(映射)与reduce(规约)

参考: monocle2

查看不同细胞群的中位基因也是一样

查看不同样品的中位基因也是一样

或者也可以



  • 单细胞数据处理小细节汇总
    答:5. 关于有些细胞属于同一个cluster但是在umap或者tsne图上相聚较远的问题:UMAP和TSNE是各自的算法在PCA降维的基础上再进行非线形降维,在二维图上把其各自算法认为相近的细胞聚在一起。但是FindClusters输入的不是UMAP或TSNE降维的数据,而是FindNeighbors的数据,而FindNeighbors输入的数据是PCA降维数据,是用另外一种算法...
  • 单细胞测序| 数据预处理常见问题
    答:Seurat3凭借其CCA降维技术,特别适合处理强相关数据,然而过度矫正可能会影响真实性的体现。如果你面对的是大样本和稀有细胞的混合,Harmony的PCA降维迭代法则更为适用,它更适用于个体差异大的临床样本,需要在矫正效果与保持原始信息之间寻找平衡。务必确保矫正的程度符合实验预期,避免矫枉过正。筛选高变基因...
  • 【单细胞测序数据分析-1】认识Seurat对象数据结构/数据格式及操作_百 ...
    答:可以通过 @active.assay 查看当前默认的assay,通过 DefaultAssay() 更改当前的默认assay。 结构 counts 存储原始数据,是稀疏矩阵 data存储logNormalize() 规范化的data。总表达式对每个单元格的要素表达式度量进行标准化,将其乘以比例因子(默认为10,000),并对结果进行对数转换 scale.data#...
  • 单细胞分析实录(6): 去除批次效应/整合数据
    答:最后还有一个有用的小技巧,在降维之后,我们可以用DimPlot()函数把结果画在二维平面上,那这些坐标信息存储在哪儿呢?存储在Embeddings(test.seu,"pca")矩阵中,harmony/tsne/umap同理。有时候我们可以把坐标和meta@data的信息提取出来,用ggplot2画图,方便修改。
  • 完整的单细胞分析流程——数据标化(normalization)
    答:假设两个细胞之间多数非DE基因之间的计数大小的任何系统性差异都代表了偏差,该偏差用于计算适当的大小因子以将其去除。 然而,由于存在大量的低计数和零计数,单细胞数据应用这些bulk归一化方法可能会有问题。为了克服这个问题,我们汇总了许多细胞的计数以进行准确的大小因子估算。然后,将基于库的大小因子“分解”为基于...
  • 完整的单细胞分析通用流程——从数据到可视化
    答:1.我们计算质量控制指标,以去除会干扰下游分析的低质量细胞。这些细胞在处理过程中可能已经损坏,或者可能没有被测序方案完全捕获。常用指标包括每个细胞的总计数,spike-in或线粒体reads的比例以及检测到的feature的数量。2.我们将计数转换为标准化的表达值,以消除特定于细胞的偏倚(例如捕获效率)。这使...
  • 单细胞测序注意事项必读
    答:13. 技术平台由于孔径限制对于细胞尺径大小有要求,一般要求在40um以下,对于肌肉细胞、脂肪细胞、神经元细胞等较大尺径细胞无法直接捕获,可以通过抽提细胞核的方法进行单细胞研究;14. 单细胞测序数据的标准分析内容只是半成品,必须经过细胞类型鉴定、细胞状态判定、轨迹推断和通路富集等一系列探索性、...
  • 单细胞RNA系列专题之一:单细胞RNA测序中质控之重要细节 (上篇)_百度...
    答:在做单细胞测序的之前,需要对细胞进行裂解。不同的细胞组织,裂解条件也会不一样。如果裂解条件过于严格,就会影响文库制备。构建文库同时加入浓度已知的spike-in,其中包括:Spike-ins 的用途 1.去除技术噪音 2.检测捕获效率 3.计算RNA的起始量 4.数据的normalization Spike-ins的问题 ...
  • 单细胞RNA系列专题之一:单细胞RNA测序中质控之重要细节 (下篇)_百度...
    答:整个单细胞分析的核心其实就是确定cell types/ lineages。而在此之前的一步就是数据质控(QC, quanlity control)。我们在得到表达矩阵之后,会做Data normalization , 基因集筛选,批次效应的去除等工作;之后用PCA, t-SNE进行降维。如果在这一过程中发现了一些问题,我们会移除掉一些细胞,然后重新质控,...
  • 6.单细胞 RNA-seq:归一化和 PCA 分析
    答:假设您正在处理 12,000 个细胞 的单细胞 RNA-seq 数据集,并且您已经量化了 20,000 个基因的表达 。计算 PC 分数后,您会看到一个 12,000 x 12,000 的矩阵,该矩阵表示有关所有细胞中相对基因表达的信息。您可以选择 PC1 和 PC2 列并以二维方式绘制它们。您还可以使用前 40 个 PC 的 PC ...