联系我们 两种不同的风光已毕harmony的多个单细胞整合

本条记会被收录于《生信手段树》公众号的《单细胞2024》专辑，况且咱们从2024启动的教程皆是基于Seurat的V5版块啦联系我们，之前依然演示了如何读取不同时势的单细胞转录组数据文献，如下所示：

初试Seurat的V5版块

使用Seurat的v5来读取多个10x的单细胞转录组矩阵

使用Seurat的v5来读取多个不是10x表率文献的单细胞面目

然则其它代码基本上就跟Seurat早期的v4莫得永诀，比如harmony整合多个单细胞样品。

但内容上Seurat有了我方的更变，官网文档：

https://satijalab.org/seurat/articles/seurat5_integration

https://satijalab.org/seurat/articles/integration_introduction

内部提到了它内置了多种整合多个单细胞样品的算法，不错 Perform streamlined (one-line) integrative analysis

Anchor-based CCA integration (method=CCAIntegration)

Anchor-based RPCA integration (method=RPCAIntegration)

Harmony (method=HarmonyIntegration)

FastMNN (method= FastMNNIntegration)

scVI (method=scVIIntegration)

上期奖号和值为97，最近十期和值分别为116 105 118 106 100 103 84 137 64 97，最近十期和值分布在64-137之间。综合分析本期预计红球和值出现在123左右。

上期龙头开出奇数号码05，近10期龙头奇数号码开出7次，偶数号码开出3次，本期优先考虑奇数号码，龙头参考05。

要是是一次性读取了多个10x样品

这个手艺，因为函数Read10X不错一次性读取多个合理的旅途，是以咱们会把多个样品就被和洽读取成为了一个寥落矩阵而不是每个样品孤苦的寥落矩阵，如下所示；

图片

和洽读取成为了一个寥落矩阵

详见：使用Seurat的v5来读取多个10x的单细胞转录组矩阵，它就不相宜走Seurat的v5的内置的多个单细胞样品的整划算法，是以咱们会先split它，代码如下所示：table(sce.all$orig.ident) 

obj = sce.all

obj[["RNA"]] <- split(obj[["RNA"]], f = obj$orig.ident)

后果如下所示，不错看到每个样品的矩阵这个手艺被上头的split函数阻隔了：

图片

split函数阻隔

接下来，如下所示走内置的harmony经过，联系我们便是IntegrateLayers函数内部的 ：obj <- NormalizeData(obj)

obj <- FindVariableFeatures(obj)

obj <- ScaleData(obj)

obj <- RunPCA(obj)

# 运行的harmony需要基于上头的PCA截至哦

sce.all <- IntegrateLayers(

object = obj, method = HarmonyIntegration,

orig.reduction = "pca", new.reduction = "harmony",

verbose = FALSE

)

sce.all@reductions$harmony

sce.all <- FindNeighbors(sce.all, reduction = "harmony", dims = 1:30)

res =  c(0.01, 0.05, 0.1, 0.2, 0.3, 0.5,0.8 )

sce.all <- FindClusters(sce.all, resolution = res )

sce.all <- RunUMAP(sce.all, reduction = "harmony", dims = 1:30)

p1 <- DimPlot( sce.all,label = T)

p1 

是很容易看到harmony被告捷运行啦。

要是是先我方跑RunHarmony函数

这个手艺又不成用split函数阻隔了的Seurat对象哦，需要最启动阿谁多个样品就被和洽读取成为了一个寥落矩阵的Seurat对象，这个手艺不错定名为 input_sce 变量，然后走底下的代码：print(dim(input_sce))

input_sce <- NormalizeData(input_sce, 

normalization.method = "LogNormalize",

scale.factor = 1e4) 

input_sce <- FindVariableFeatures(input_sce)

input_sce <- ScaleData(input_sce)

input_sce <- RunPCA(input_sce, features = VariableFeatures(object = input_sce))

seuratObj <- RunHarmony(input_sce, "orig.ident")

names(seuratObj@reductions)

seuratObj <- RunUMAP(seuratObj,  dims = 1:15, 

reduction = "harmony")

# p = DimPlot(seuratObj,reduction = "umap",label=T ) 

# ggsave(filename='umap-by-orig.ident-after-harmony',plot = p)

input_sce=seuratObj

input_sce <- FindNeighbors(input_sce, reduction = "harmony",

dims = 1:15)  

这个后果跟前边的IntegrateLayers函数是一趟事，要是你这个手艺的Seurat对象内部的矩阵是按照样品阻隔的， 就需要先统一智商走RunHarmony函数，它来自于harmony这个r包啦。

常常是不淡薄走内置的harmony经过

其实我就不应该先容这个IntegrateLayers函数的，因为它需要split阿谁矩阵，这么的话后头的好多分析皆会有问题，比如咱们跑 cosg 函数针对阿谁矩阵去找marker就会遭遇报错：#  remotes::install_github('genecell/COSGR')

#  genexcell <- Seurat::GetAssayData(object = object[[assay]],slot = slot)

marker_cosg <- cosg(

sce.all.int,

groups='all',

assay='RNA',

slot='data',

mu=1,

n_genes_user=100)

如下所示：Error in `Seurat::GetAssayData()`:

! GetAssayData doesn't work for multiple layers in v5 assay.

Run `rlang::last_trace()` to see where the error occurred.

> rlang::last_trace()

<error  you can run 'object 

Error in `Seurat::GetAssayData()`:

! GetAssayData doesn't work for multiple layers in v5 assay.

---

Backtrace:

▆

1. └─COSG::cosg(...)

2.   ├─Seurat::GetAssayData(object = object[[assay]], slot = slot)

3.   └─SeuratObject:::GetAssayData.StdAssay(object = object[[assay]], slot = slot)

Run rlang::last_trace(drop = FALSE) to see 1 hidden frame.

> 

要是要科罚这个报错，就需要最初把前边的split的矩阵重新joint且归，又是艰巨的事情！！！

文末友情宣传

热烈淡薄你保举给身边的博士后以及年青生物学PI，多小数数据理解，让他们的科研上一个台阶：

生物信息学马拉松讲课（买一得五） ，你的生物信息学初学课，春节前临了一次

时隔5年，咱们的生信手段树VIP学徒持续招生啦

144线程640Gb内存作事器分享一年仍然是仅需800联系我们

• 联系我们
• 不同
• 风光
• 两种
• harmony