查看原文
其他

我确实不知道harmony是全球首款性爱机器人

生信技能树 生信技能树 2022-08-15

多个病人的多个样本进行单细胞转录组测序是非常正常的,比如下面的:

  • 发表在 Nat Med. 2018 Aug; 题目是:Phenotype molding of stromal cells in the lung tumor microenvironment. 共选取5例病人的共19个样本,通过10×genomics单细胞转录组测序探索基质细胞的亚群分类、基因功能(信号通路)、关键marker基因和临床预后,共鉴定出52个基质细胞亚群,
  • 发表在  Nature Medicine  (2018) ,标题是:Single-cell profiling of breast cancer T cells reveals a tissue-resident memory subset associated with improved prognosis  作者从3个乳腺癌患者体内通过FACS筛选到乳腺癌中肿瘤浸润淋巴细胞(TILs),使用商业仪器10X来做单细胞转录组,得到6,311个T细胞数据。
  • 还有我们在单细胞天地分享的scRNA-seq揭示胰腺导管腺癌的瘤内异质性和恶性进展 ,有24个原发性PDAC肿瘤病人样本及11个对照胰腺样本(3例非胰腺肿瘤患者和8例非恶性胰腺肿瘤患者样本)

这些2018年的研究,大多数都是百万经费起的项目,所以基本上是CNS级别的成果。现在想发普通的单细胞文章,也是得做多个样本了,就面临如何整合的问题,其中最出名的当然是Seurat包的CCA方法了,具体多火呢,发了才一年,引用就快破千!

这样的单细胞转录组测序整合算法超级多,我想你应该是不会一一细看的,我就列出来而已:

  • MNNcorrect (https://doi.org/10.1038/nbt.4091)
  • CCA + anchors (Seurat v3) (https://doi.org/10.1101/460147)
  • CCA + dynamic time warping (Seurat v2) (https://doi.org/10.1038/nbt.4096)
  • LIGER (https://doi.org/10.1101/459891)
  • Harmony (https://doi.org/10.1101/461954)
  • Conos(https://doi.org/10.1101/460246)
  • Scanorama(https://doi.org/10.1101/371179)
  • scMerge(https://doi.org/10.1073/pnas.1820006116)

今天就从简单的Seurat包的merge函数到,CCA + anchors (Seurat v3),再到Harmony 给大家介绍一下哈。

因为一直听说harmony效果更佳,今天才有空所以想体验一下,比如 2019的 Cell的文章《Landscap and Dynamics of Single Immune Cells in Hepatocellular Carcinoma》, 整合两类数据(包括SMART-seq2和10X),就使用的是 harmony算法。

结果在“百度百科”查到它:据英国《每日邮报》2017年4月4日报道,RealDoll公司设计出世界上首款性爱机器人,可以和用户产生情感交流。Harmony 2.0配有苏格兰口音,即将上市,售价预计为1万美元。(如果有关于Harmony的用户体验欢迎交流哈!)

先看看简单的Seurat包的merge函数

首先保证你的 folder 变量存储着全部的10x单细胞转录组数据的cellranger结果哈;

# folder 这个变量,存储着多个10x单细胞转录组数据的cellranger结果
library(Seurat)
sceList = lapply(folders,function(folder){ 
  CreateSeuratObject(counts = Read10X(folder), 
                               project = folder )
})
 
sce.big <- merge(sceList[[1]], 
                 y = c(sceList[[2]],sceList[[3]],sceList[[4]],sceList[[5]],
                       sceList[[6]],sceList[[7]],sceList[[8]]), 
                 add.cell.ids = folders, 
                 project = "mouse8")
sce.big
table(sce.big$orig.ident)
save(sce.big,file = 'sce.big.merge.mouse8.Rdata')

然后Seurat包的CCA整合

# folder 这个变量,存储着多个10x单细胞转录组数据的cellranger结果
library(Seurat)
sceList = lapply(folders,function(folder){ 
  CreateSeuratObject(counts = Read10X(folder), 
                               project = folder )
})
sceList
# 标注化和寻找 VariableFeatures
for (i in 1:length(sceList)) {
  print( Sys.time() )
  sceList[[i]] <- NormalizeData(sceList[[i]], verbose = FALSE)
  sceList[[i]] <- FindVariableFeatures(sceList[[i]], selection.method = "vst"
                                       nfeatures = 2000, verbose = FALSE)
}

sceList
Sys.time()
#  整合 二个样本
sce.anchors <- FindIntegrationAnchors(object.list = sceList, dims = 1:30)
# load('sce.anchors.RDATA')
Sys.time()
sce.integrated <- IntegrateData(anchorset = sce.anchors, dims = 1:30)
Sys.time()



library(ggplot2)
library(cowplot)
# switch to integrated assay. The variable features of this assay are automatically
# set during IntegrateData
DefaultAssay(sce.integrated) <- "integrated"
# 接下来的分析,都基于 sce.integrated 啦

接着是MNNs算法

Seurat主要是处理10x单细胞转录组数据,而10x仪器商业上的成功可以说是成就了Seurat包,另外一个比较火的多个样本单细胞转录组数据整合算法是mutual nearest neighbors (MNNs)

因为并不是必须,所以我一直没有使用它!单细胞转录组数据的整合分析才是主流,所以我在生信技能树的多个单细胞转录组样本的数据整合之CCA-Seurat包 是非常值得大家细看的,而且使用scran包的MNN算法来去除多个单细胞转录组数据批次效应 也展现了一个很好的例子。

最后是我们今天的主角harmony

发表这个harmony算法的文章是2019年11月的《Fast, sensitive and accurate integration of single-cell data with Harmony》,目前仍然是一个GitHub包,在 https://github.com/immunogenomics/harmony 所以安装它的方式是:

library(devtools)
install_github("immunogenomics/harmony")

前面的Seurat标准流程仍然是一样的:

sce <- NormalizeData(sce, normalization.method =  "LogNormalize",  
                     scale.factor = 1e4)
GetAssay(sce,assay = "RNA")
sce <- FindVariableFeatures(sce, 
                            selection.method = "vst", nfeatures = 2000)  
sce <- ScaleData(sce) 
sce <- RunPCA(object = sce, pc.genes = VariableFeatures(sce)) 
DimHeatmap(sce, dims = 1:12, cells = 100, balanced = TRUE)
ElbowPlot(sce) 

这个时候的sce对象,里面存储这多个10x样品的单细胞结果,以 orig.ident 列作为区分。

接下来开始引入harmony步骤,超级简单:

library(harmony)
seuratObj <- RunHarmony(sce, "orig.ident")
names(seuratObj@reductions)
seuratObj <- RunUMAP(seuratObj,  dims = 1:15
                     reduction = "harmony")
DimPlot(seuratObj,reduction = "umap",label=T ) 

对harmony整合好的对象进行聚类分群:

sce=seuratObj
sce <- FindNeighbors(sce, reduction = "harmony",dims = 1:15)
sce <- FindClusters(sce, resolution = 0.1)
table(sce@meta.data$seurat_clusters)
DimPlot(sce,reduction = "umap",label=T)  
DimPlot(sce,reduction = "umap",label=T,
        group.by = 'orig.ident'

可以看到如下所示的结果,上图是整合前,下图整合后 :

上图是整合前,下图整合后

如果你感兴趣这个harmony技术的细节,建议看前两天在《单细胞天地》的 《单细胞分析十八般武艺1:harmony》

类似的容易让人想歪的生物信息学词条好像还不少

比如:julia 这个编程语言,如果你 直接搜索它!

julia这个女明星

哈哈,群体的智慧是无穷的,欢迎各路神通广大的朋友留言给出类似的例子!

您可能也对以下帖子感兴趣

文章有问题?点此查看未经处理的缓存