标准流程参数调整 dims

2023-05-14 本文已影响0人 oceanandshore

seurat去批次标准流程

https://www.jianshu.com/p/e1dbe1062419
https://www.jianshu.com/p/4e748d3715a8

library(harmony)
library(devtools)
library(Seurat)
library(tidyverse)
library(patchwork)
library(ggplot2)
library(ggraph)
library(sctransform)
library(glmGamPoi)

rm(list=ls())

###1、样本读取，构建对象


Use.data <- Read10X(data.dir = "D:/atrial fibrillation 2022.11.03/GSE1109_RAW1/Use")
Use <- CreateSeuratObject(counts = Use.data, min.cells=3,
                          min.features = 200,project = "Use")
Use

Withdrawal.data <- Read10X(data.dir = "D:/atrial fibrillation 2022.11.03/GSE1109_RAW1/Withdrawal")
Withdrawal <- CreateSeuratObject(counts = Withdrawal.data, min.cells=3,
                                 min.features = 200,project = "Withdrawal")
Withdrawal


###2、分组   https://zhuanlan.zhihu.com/p/465893822


Use$group<-"Use"
Withdrawal$group<-"Withdrawal"


###3、设置dataset_list

dataset_list <- list(Use,Withdrawal)


###4、Qc+标准化

dataset_list <- lapply(dataset_list,function(object) {
  object <- PercentageFeatureSet(object, pattern = "^MT-", col.name = "percent.mt")
  object <- subset(object,subset = nFeature_RNA > 800 & nFeature_RNA < 2800 &nCount_RNA >1000 &nCount_RNA <20000 & percent.mt < 13)#参数过滤细胞
  object <- NormalizeData(object, verbose = FALSE) 
  object <- FindVariableFeatures(object, selection.method = "vst", 
                                 nfeatures = 2000, verbose = FALSE)
})

table(integrated@meta.data$orig.ident) 

###5、整合数据去批次IntegrateData

integration_features <- SelectIntegrationFeatures(object.list = dataset_list)
integration_anchors <- FindIntegrationAnchors(object.list = dataset_list, anchor.features = integration_features)
integrated <- IntegrateData(anchorset = integration_anchors)


###6、switch to integrated assay. 
DefaultAssay(integrated) <- "integrated"


###7、后续
integrated <- ScaleData(integrated, verbose = FALSE)
saveRDS(integrated, file = "./testtesttesttesttesttesttesttesttesttest.rds")


##########################################################################

rm(list=ls())
integrated = readRDS( file = "./testtesttesttesttesttesttesttesttesttest.rds")

integrated <- RunPCA(integrated, npcs = 50, verbose = FALSE)
integrated <- RunUMAP(integrated, reduction = "pca", dims = 1:50)
integrated <- FindNeighbors(integrated, reduction = "pca", dims = 1:50)
integrated <- FindClusters(integrated, resolution = 0.5)


plot4 = DimPlot(integrated, reduction = "umap", group.by='orig.ident') 
plot6<- DimPlot(integrated, reduction = "umap", label = TRUE, repel = TRUE)
#combinate
plotc <- plot4+plot6
plotc

library(Seurat)
library(Azimuth)
library(SeuratData)
library(patchwork)

# Install the PBMC systematic comparative analyis (pmbcsca) dataset
InstallData("pbmcsca")

# returns a Seurat object named pbmcsca
pbmcsca <- LoadData("pbmcsca")

# The RunAzimuth function can take a Seurat object as input
integrated <- RunAzimuth(integrated, reference = "pbmcref")


p1 <- DimPlot(integrated, group.by = "predicted.celltype.l2", label = TRUE, label.size = 3) + NoLegend()
p3 <- DimPlot(integrated,label = TRUE,repel = TRUE, label.size = 3) 
p1+p3

dim=20 resolution = 0.5 第7群分不开，AZI注释显示第7群是dnT细胞，查了前几个特异性marker也确实是dnT细胞。但是明显dnT细胞的比例没这么高，结合第7群中间不是那么紧密，应该是第7群没分开的原因。