如果想画空间转录组的细胞之间的相关性，需要先使用cell2location计算细胞丰度
单细胞、空转联合分析| 揭示结直肠癌中FAP+成纤维细胞和SPP1+巨噬细胞的相互作用
加载数据 #request 2.libPaths(c( "/home/data/t040413/R/x86_64-pc-linux-gnu-library/4.2","/home/data/t040413/R/yll/usr/local/lib/R/site-library","/home/data/refdir/Rlib/", "/usr/local/lib/R/library"))setwd('/home/data/t040413/silicosis/spatial/')dir.create('/home/data/t040413/silicosis/spatial/celltype_inspot_enrichments')setwd('/home/data/t040413/silicosis/spatial/celltype_inspot_enrichments')getwd()library(Seurat)library(dplyr)#d.all=readRDS("~/silicosis/spatial/integrated_slides/integrated_slides.rds")load("/home/data/t040413/silicosis/spatial_transcriptomics/silicosis_ST_harmony_SCT_r0.5.rds")dim(d.all)SpatialFeaturePlot(d.all,features = c('Grem1','Inmt','Runx1'))d.all$stim=factor(d.all$stim,levels = c( "NS_7","SiO2_7","NS_56","SiO2_56"))Idents(d.all)= factor( Idents(d.all),levels = c( "NS_7","SiO2_7","NS_56","SiO2_56"))#save(d.all,file ="/home/data/t040413/silicosis/spatial_transcriptomics/silicosis_ST_harmony_SCT_r0.5.rds" )myimages=d.all@images# Reordering using indexingnew_order <- c("NS_7", "SiO2_7", "NS_56","SiO2_56")myimages <- myimages[new_order]d.all@images=myimages
DefaultAssay(d.all)DefaultAssay(d.all)="Spatial"#visium_slides=SplitObject(object = d.all,split.by = "stim")names(d.all);dim(d.all)d.all@meta.data %>%head()head(colnames(d.all))#给d.all 添加meta信息------adata_obs=read.csv("~/silicosis/spatial/adata_obs.csv")head(adata_obs)mymeta= paste0(d.all@meta.data$orig.ident,"_",colnames(d.all)) %>% gsub("-.*","",.) # %>% head()head(mymeta)tail(mymeta)#掉-及其之后内容adata_obs$col= adata_obs$spot_id %>% gsub("-.*","",.) # %>% head()head(adata_obs)rownames(adata_obs)=adata_obs$coladata_obs=adata_obs[mymeta,]head(adata_obs)identical(mymeta,adata_obs$col)d.all=AddMetaData(d.all,metadata = adata_obs)DefaultAssay(d.all)='SCT'DimPlot(d.all,group.by = 'stim')DimPlot(d.all,label = TRUE)dim(d.all)dim(d.all@assays$SCT)dim(d.all@assays$Spatial)##2 添加细胞丰度信息-------#major_abundance1=read.csv("~/silicosis/spatial_c2l_results/sc_silicosis_20231006/cell2location_map/q05_cell_abundance_w_sf.csv")
加载cell2location的细胞丰度结果
##2 添加细胞丰度信息-------major_abundance1=read.csv("~/silicosis/spatial_c2l_results/silicosis_mapped_after_kaiti/sc_silicosis_analysis/cell2location_map/q05_cell_abundance_w_sf.csv")head(major_abundance1)rownames(major_abundance1)= major_abundance1$spot_id %>%stringr::str_split(pattern = "-",simplify = T) %>%.[,1]major_abundance=major_abundance1[,-1]head(major_abundance)#3对齐细胞丰度文件与seurat的barcodeshead(colnames(d.all))head(d.all@meta.data)major_abundance=major_abundance[d.all$col,]head(major_abundance)#把细胞丰度信息添加到seurat对象中d.all=AddMetaData(d.all,metadata = major_abundance)SpatialFeaturePlot(d.all,features = 'Inmt.fibroblast')SpatialFeaturePlot(d.all,features = 'Grem1.fibroblast')exprSet=as.matrix(major_abundance)head(exprSet)
准备画图数据
selected_Genes=colnames(exprSet)FUN = function(selected_Genes,exprSet){if (T) {y=as.numeric(exprSet[,selected_Genes])colnames <- colnames(exprSet)cor_data_df <- data.frame(colnames)head(cor_data_df)for (i in 1:length(colnames)){test <- cor.test(as.numeric(exprSet[,i]),y,type="spearman")cor_data_df[i,2] <- test$estimatecor_data_df[i,3] <- test$p.valuecor_data_df[i,4]=selected_Genes}head(cor_data_df)names(cor_data_df) <- c("symbol","spearman_correlation","spearman_pvalue",'symbol2')return(cor_data_df)#head(cor_data_df)}}cor_data_list=lapply(selected_Genes,FUN ,exprSet = exprSet)cor_data_listcor_data_df=do.call(rbind,cor_data_list)head(cor_data_df)write.csv(cor_data_df,file = "./cor_celltype_inspot_data_df.csv")sio2_7=exprSet[stringr::str_detect(string = rownames(exprSet),pattern = 'SiO2_7') ,]sio2_56=exprSet[stringr::str_detect(string = rownames(exprSet),pattern = 'SiO2_56') ,]sio2=exprSet[stringr::str_detect(string = rownames(exprSet),pattern = 'SiO2') ,]ns=exprSet[stringr::str_detect(string = rownames(exprSet),pattern = 'NS') ,]ns_7=exprSet[stringr::str_detect(string = rownames(exprSet),pattern = 'NS_7') ,]ns_56=exprSet[stringr::str_detect(string = rownames(exprSet),pattern = 'NS_56') ,]head(exprSet)print(getwd())#install.packages("dplyr")#install.packages("ggplot2")library(ggstatsplot)
画图
p1=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(exprSet),x = "Grem1.fibroblast",y='AM3' ) ;p1p2=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(exprSet),x = 'AM3',y= "Grem1.fibroblast") ;p2getwd()ggsave(filename = "./celltype_inspot_correlation_grem1_am3-.png",plot = p2,width = 16,height = 8,dpi = 900,limitsize = FALSE)library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_grem1_am3.pdf",plot = p1,width = 9,height = 9)ggsave(filename = "./celltype_inspot_correlation_grem1_am32.pdf",plot = p1,width = 9,height = 8)p2=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(exprSet),x = 'Inmt.fibroblast',y= "Grem1.fibroblast") ;p2p3=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(exprSet),x = 'Inmt.fibroblast',y= "AM3") ;p3#p=VlnPlot(sce.all.int,group.by = "group")#sio2_7----head(exprSet)p2=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(sio2_7),x = "Grem1.fibroblast",y='AM3');p2library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_grem1_am3—sio2_7.pdf",plot = p2,width = 9,height = 9)p3=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(sio2_56),x = "Grem1.fibroblast",y='AM3');p3library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_grem1_am3—sio2_56.pdf",plot = p3,width = 9,height = 9)p4=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(sio2),x = "Grem1.fibroblast",y='AM3');p4library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_grem1_am3—sio2.pdf",plot = p4,width = 9,height = 9)p5=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(ns),x = "Grem1.fibroblast",y='AM3');p5library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_grem1_am3—ns.pdf",plot = p5,width = 9,height = 9)p6=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(ns_7),x = "Grem1.fibroblast",y='AM3')library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_grem1_am3—ns7.pdf",plot = p6,width = 9,height = 9)p7=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(ns_56),x = "Grem1.fibroblast",y='AM3')library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_grem1_am3—ns56.pdf",plot = p7,width = 9,height = 9)head(exprSet)p8=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(ns_56),x = "AM3",y='Inmt.fibroblast')library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_INMT_am3—ns56.pdf",plot = p8,width = 9,height = 9)p9=ggstatsplot::ggscatterstats(data = as.data.frame.matrix(sio2_56),x = "AM3",y='Inmt.fibroblast')library(ggplot2)ggsave(filename = "./celltype_inspot_correlation_INMT_am3—sio2_56.pdf",plot = p9,width = 9,height = 9)table(exprSet[,"Grem1.fibroblast"]==0)table(exprSet[,'AM3']==0)head(exprSet)library(dplyr)long_exprset=exprSet %>% as.data.frame() %>% mutate("barcodes"=rownames(.)) %>%tidyr:: pivot_longer(cols = -barcodes,# cols = colnames(exprSet),names_to = "cell types",values_to = "cell abundance")long_exprset2=long_exprset[long_exprset$`cell abundance`>0.2,]long_exprset2wide_exprset <- long_exprset2 %>%tidyr:: pivot_wider(names_from = "cell types",values_from = "cell abundance")wide_exprsetggstatsplot::ggscatterstats(data = as.data.frame.matrix(wide_exprset),x = "Grem1.fibroblast",y='AM3')ggplot2::ggplot(data = long_exprset,mapping = )table(long_exprset$`cell types`)# 创建散点图scatter_plot <- ggplot(long_exprset$`cell types` %in% c(""), aes(x = `cell abundance`, y = `cell types`, color = `cell types`)) +geom_point() +labs(x = "cell abundance", y = "cell types") +theme_minimal() +geom_smooth(method = "lm", se = FALSE)# 添加相关性曲线
