作者，追风少年i
新的一周，新的开始，接上一篇，10X空间转录组数据分析之细胞的空间依赖性，这一篇我们来分析空间细胞类型和空间通路分布的依赖性，如下图：
图片.png
我们的目的就是研究空间细胞类型分布是否与生物学通路的分布存在强相关性。
并且我们要用代码来实现一下
library(tidyverse)
library(Seurat)
library(mistyR)
source("./misty_utilities.R")
source的misty_utilities.R内容在文章10X空间转录组数据分析之细胞的空间依赖性最下面
我们需要准备的内容
- 单细胞空间联合分析的矩阵（Seurat、cell2location均可）
- 空间每个spot的通路分析结果（富集分析，或者其他方式的打分矩阵）
定义共定位函数
# Pipeline definition:
run_colocalization <- function(slide,
assay,
useful_features,
useful_features_ct,
out_label,
misty_out_alias = "./results/tissue_structure/misty/pathway_map/pm_") {
# Define assay of each view ---------------
view_assays <- list("main" = assay,
"juxta" = assay,
"para" = assay,
"intra_ct" = "c2l",
"para_ct" = "c2l")
# Define features of each view ------------
view_features <- list("main" = useful_features,
"juxta" = useful_features,
"para" = useful_features,
"intra_ct" = useful_features_ct,
"para_ct" = useful_features_ct)
# Define spatial context of each view -----
view_types <- list("main" = "intra",
"juxta" = "juxta",
"para" = "para",
"intra_ct" = "intra",
"para_ct" = "para")
# Define additional parameters (l in case of paraview,
# n of neighbors in case of juxta) --------
view_params <- list("main" = NULL,
"juxta" = 5,
"para" = 15,
"intra_ct" = NULL,
"para_ct" = 15)
misty_out <- paste0(misty_out_alias,
out_label, "_", assay)
run_misty_seurat(visium.slide = slide,
view.assays = view_assays,
view.features = view_features,
view.types = view_types,
view.params = view_params,
spot.ids = NULL,
out.alias = misty_out)
return(misty_out)
}
读取空间的rds文件和spot打分注释文件
slide = readRDS(spatialrds)
en_anno = read.csv(enrichment_spatial.csv) ###这里大家需要自己分析
这里的通路富集分析推荐大家使用R包progeny，PROGENy is resource that leverages a large compendium of publicly available signaling perturbation experiments to yield a common core of pathway responsive genes for human and mouse. These, coupled with any statistical method, can be used to infer pathway activities from bulk or single-cell transcriptomics.
library(progeny)
if(species == "mouse"){
model <- progeny::getModel(organism = "Mouse", top = top)
}else if(species == "human"){
model <- progeny::getModel(organism = "Human", top = top)
富集分析推荐使用AUC或者GSVA，这里就不做了，直接拿着分析结果往下
progeny_scores <- scale(t(progeny_scores)) ####分数的标准化操作
slide[['progeny']] <- CreateAssayObject(counts = t(progeny_scores))
其二是单细胞空间联合分析的矩阵
anno = read.csv(sp_sc.csv,header = T,row.names = 1)
slide[['c2l']] <- CreateAssayObject(counts = t(anno))
数据分析
assay_label <- "progeny"
assay <- assay_label
DefaultAssay(slide) <- assay
useful_features <- rownames(slide)
useful_features <- useful_features[! useful_features %in% c("TNFa")] ###通路特征
useful_features_ct <- rownames(GetAssayData(slide, assay = "c2l")) ###细胞类型特征
useful_features_ct <- useful_features_ct[! useful_features_ct %in% "prolif"]
mout <- run_colocalization(slide = slide,
useful_features = useful_features,
useful_features_ct = useful_features_ct,
out_label = slide_id,
assay = assay,
misty_out_alias = "./results/tissue_structure/misty/pathway_map/pm_")
misty_res_slide <- collect_results(mout)
数据分析结束后进行可视化
plot_folder <- paste0(mout, "/plots")
system(paste0("mkdir ", plot_folder))
pdf(file = paste0(plot_folder, "/", slide_id, "_", "summary_plots.pdf"))
mistyR::plot_improvement_stats(misty_res_slide)
mistyR::plot_view_contributions(misty_res_slide)
mistyR::plot_interaction_heatmap(misty_res_slide, "intra", cutoff = 0)
mistyR::plot_interaction_communities(misty_res_slide, "intra", cutoff = 0.5)
mistyR::plot_interaction_heatmap(misty_res_slide, "juxta_5", cutoff = 0)
mistyR::plot_interaction_communities(misty_res_slide, "juxta_5", cutoff = 0.5)
mistyR::plot_interaction_heatmap(misty_res_slide, "para_15", cutoff = 0)
mistyR::plot_interaction_communities(misty_res_slide, "para_15", cutoff = 0.5)
mistyR::plot_interaction_heatmap(misty_res_slide, "intra_ct", cutoff = 0)
mistyR::plot_interaction_heatmap(misty_res_slide, "para_ct_15", cutoff = 0)
dev.off()
最后分析得到细胞类型分布与空间分布的依赖性热图
图片.png
好了，已经分享给大家了，生活很好，有你更好，百度文库出现了大量抄袭我的文章，对此我深表无奈，我写的文章，别人挂上去赚钱，抄袭可耻，挂到百度文库的人更可耻
