From b88b66137826550539da66a774a9b87329e8e647 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=E7=94=9F=E4=BF=A1=E5=88=86=E6=9E=90?= Date: Wed, 11 Dec 2024 02:53:16 +0800 Subject: [PATCH] =?UTF-8?q?=E6=B7=BB=E5=8A=A0=2010xscRNA-seq.R?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: 生信分析 --- 10xscRNA-seq.R | 114 +++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 114 insertions(+) create mode 100644 10xscRNA-seq.R diff --git a/10xscRNA-seq.R b/10xscRNA-seq.R new file mode 100644 index 0000000..49019c7 --- /dev/null +++ b/10xscRNA-seq.R @@ -0,0 +1,114 @@ +## Written By LiShang +## Notice: 仅适用于10x平台 + +#加载依赖包 +library(Seurat) +library(sctransform) +library(dplyr) +library(ggplot2) + +########################读取样本并生成Seurat Object########################## + +#获取所有样本的存放目录 +setwd("~/scRNA-seq") +sample_dirs <- list.dirs(full.names = FALSE)[-1] + +#同时读取所有样本, 并转化成Seurat Object +sample_list <- mapply(function(dir) { + Read10X(data.dir = dir) %>% + CreateSeuratObject(project = dir, min.cells = 3, min.features = 200) +}, dir = sample_dirs, SIMPLIFY = FALSE) +#合并所有样本 +merged <- merge(sample_list[[1]], sample_list[-1], add.cell.ids = sample_dirs) + +#提取样本的分组信息, 本模板分为cancer和normal两个组 +merged[["sample.grp"]] <- merged@meta.data$orig.ident %>% + strsplit(split = "_") %>% + sapply(function(x) x[2]) %>% + as.factor() + +#################################质量控制#################################### + +#计算线粒体转录本占比 +merged[["percent.mt"]] <- PercentageFeatureSet(merged, pattern = "^MT-") +#细胞过滤 +merged <- merged %>% + subset(subset = nFeature_RNA > 200 & nFeature_RNA < 2500) %>% + subset(subset = percent.mt < 5) + +###############################标准分析流程################################## + +#SCTransform标准化、线性降维 +merged <- merged %>% + SCTransform(vars.to.regress = "percent.mt", verbose = FALSE) %>% + RunPCA(verbose = FALSE) + +#多样本数据整合, harmony去批次 +merged <- IntegrateLayers( + object = merged, + method = HarmonyIntegration, + normalization.method = "SCT", + orig.reduction = "pca", + new.reduction = "harmony", + verbose = FALSE +) +merged[["RNA"]] <- JoinLayers(merged[["RNA"]]) + +#聚类、非线性降维 +merged <- merged %>% + FindNeighbors(dims = 1:30, reduction = "harmony", verbose = FALSE) %>% + FindClusters(resolution = 2, verbose = FALSE) %>% + RunUMAP(dims = 1:30, reduction = "harmony", verbose = FALSE) + +###############################细胞类型注释################################## + +DimPlot(merged, reduction = "umap", label = TRUE) +DotPlot(merged, c("EPCAM", "CDH1", "GATA3", "KRT18", #Epithelial + "ESR1", "PGR", "ERBB2", #Breast Cancer + "PDGFRA", "COL1A1", "FAP", "DCN", "VIM", #Fibroblast + "CD34", "VWF", "CDH5", "ENG", "PECAM1", #Endothelial + "PTPRC", #CD45: Immune + "CD2", "CD3D", "CD3E", "CD3G", #T cell + "MS4A1", "CD79A", "CD79B", #B cell + "KLRF1", "KLRD1", "NKG7", "GNLY" #NK cell + ), cols = c("yellow","blue") +) + theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1)) + +new_cluster_ids <- c( + "T cell", + "malignant", + "T cell", + "T cell", + "malignant", + "malignant", + "malignant", + "malignant", + "fibroblast", + "malignant", + "fibroblast", + "malignant", + "malignant", + "malignant", + "malignant", + "endothelial", + "fibroblast", + "T cell", + "fibroblast", + "malignant", + "malignant", + "malignant", + "malignant", + "B cell" +) +names(new_cluster_ids) <- levels(merged) +merged <- RenameIdents(merged, new_cluster_ids) + +###############################数据分析作图################################## + +merged <- PrepSCTFindMarkers(merged) +markers <- FindMarkers(merged, ident.1 = "malignant") + +FeaturePlot(merged, features = c("CREB1","CREB3")) +DotPlot(merged, features = c("CREB1","CREB3")) +VlnPlot(merged, features = c("CREB1")) +