jamesbang

V1

2022/10/23阅读：72主题：雁栖湖

# 🤗 Seurat | 超好用的单细胞测序数据合并（3'和5'数据合并）（一）

## 1写在前面

• `Harmony`;
• `rliger`;
• `Seurat`

• `3'``5'`不同`datasets`的合并；
• 整合只有部分重叠的`datasets`，（举个栗子🌰：全血`scRNAseq`数据和`3'PBMC`数据的合并。

🌟 在`Seurat`包中提供了一种叫`canonical correlation analysis` (`CCA`)的方法进行合并。

## 2用到的包

``rm(list = ls())library(Seurat)library(SeuratDisk)library(SeuratWrappers)library(patchwork)library(harmony)library(rliger)library(RColorBrewer)library(tidyverse)library(reshape2)library(ggsci)library(ggstatsplot)``

## 3示例数据

``matrix_3p <- Read10X_h5("./3p_pbmc10k_filt.h5",use.names = T)matrix_5p <- Read10X_h5("./5p_pbmc10k_filt.h5",use.names = T)\$`Gene Expression`srat_3p   <- CreateSeuratObject(matrix_3p,project = "pbmc10k_3p")srat_5p   <- CreateSeuratObject(matrix_5p,project = "pbmc10k_5p")srat_3psrat_5p``

Note! `5' datset`中还有一个`assay`，即`VDJ data`。🤜

## 4线粒体与核糖体基因的计算

``srat_3p[["percent.mt"]]  <- PercentageFeatureSet(srat_3p, pattern = "^MT-")srat_3p[["percent.rbp"]] <- PercentageFeatureSet(srat_3p, pattern = "^RP[SL]")srat_5p[["percent.mt"]]  <- PercentageFeatureSet(srat_5p, pattern = "^MT-")srat_5p[["percent.rbp"]] <- PercentageFeatureSet(srat_5p, pattern = "^RP[SL]")``

3' dataset可视化 🐶

``VlnPlot(srat_3p, ncol = 4,        features = c("nFeature_RNA","nCount_RNA","percent.mt","percent.rbp"))``

5' dataset可视化 🐶

``VlnPlot(srat_5p, features = c("nFeature_RNA","nCount_RNA","percent.mt","percent.rbp"), ncol = 4)``

``table(rownames(srat_3p) %in% rownames(srat_5p)) ``

## 5过滤

``srat_3p <- subset(srat_3p, subset = nFeature_RNA > 500 & nFeature_RNA < 5000 & percent.mt < 15)srat_5p <- subset(srat_5p, subset = nFeature_RNA > 500 & nFeature_RNA < 5000 & percent.mt < 10)``

## 6合并

### 6.1 转为list

``pbmc_list <- list()pbmc_list[["pbmc10k_3p"]] <- srat_3ppbmc_list[["pbmc10k_5p"]] <- srat_5pfor (i in 1:length(pbmc_list)) {  pbmc_list[[i]] <- NormalizeData(pbmc_list[[i]], verbose = F)  pbmc_list[[i]] <- FindVariableFeatures(pbmc_list[[i]], selection.method = "vst", nfeatures = 2000, verbose = F)}``

### 6.2 寻找合并Anchors

``pbmc_anchors <- FindIntegrationAnchors(object.list = pbmc_list, dims = 1:30)``

### 6.3 合并

``pbmc_seurat <- IntegrateData(anchorset = pbmc_anchors, dims = 1:30)rm(pbmc_list)rm(pbmc_anchors)``

## 7合并前后的比较

### 7.1 查看信息

``pbmc_seurat``

### 7.2 合并前

``DefaultAssay(pbmc_seurat) <- "RNA"pbmc_seurat <- NormalizeData(pbmc_seurat, verbose = F)pbmc_seurat <- FindVariableFeatures(pbmc_seurat, selection.method = "vst", nfeatures = 2000, verbose = F)pbmc_seurat <- ScaleData(pbmc_seurat, verbose = F)pbmc_seurat <- RunPCA(pbmc_seurat, npcs = 30, verbose = F)pbmc_seurat <- RunUMAP(pbmc_seurat, reduction = "pca", dims = 1:30, verbose = F)### 可视化p1 <- DimPlot(pbmc_seurat,reduction = "umap") +   scale_color_npg()+  plot_annotation(title = "10k 3' PBMC and 10k 5' PBMC cells, before integration")p2 <- DimPlot(pbmc_seurat, reduction = "umap", split.by = "orig.ident") +  scale_color_npg()+  NoLegend()p1 + p2``

### 7.3 合并后

``DefaultAssay(pbmc_seurat) <- "integrated"pbmc_seurat <- ScaleData(pbmc_seurat, verbose = F)pbmc_seurat <- RunPCA(pbmc_seurat, npcs = 30, verbose = F)pbmc_seurat <- RunUMAP(pbmc_seurat, reduction = "pca", dims = 1:30, verbose = F)### 可视化p1 <- DimPlot(pbmc_seurat, reduction = "umap") +   scale_color_npg()+  plot_annotation(title = "10k 3' PBMC and 10k 5' PBMC cells, after integration")p2 <- DimPlot(pbmc_seurat, reduction = "umap", split.by = "orig.ident") +  scale_color_npg()+  NoLegend()p1 + p2``

## 8聚类

### 8.1 寻找clusters

``pbmc_seurat <- FindNeighbors(pbmc_seurat, dims = 1:30, k.param = 10, verbose = F)pbmc_seurat <- FindClusters(pbmc_seurat, verbose = F)ncluster <- length(unique(pbmc_seurat[[]]\$seurat_clusters))mycol <- colorRampPalette(brewer.pal(8, "Set2"))(ncluster)DimPlot(pbmc_seurat,        cols = mycol,        label = T,repel = T) +   NoLegend()``

### 8.2 具体查看及可视化

``count_table <- table(pbmc_seurat@meta.data\$seurat_clusters, pbmc_seurat@meta.data\$orig.ident)count_table#### 可视化count_table %>%   as.data.frame() %>%   ggbarstats(x = Var2,              y = Var1,             counts = Freq)+  scale_fill_npg()``

##### jamesbang
V1

wx🔍: Grassssss 卷起来了