方法更新

admin管理员组

文章数量:1440027

方法更新

作者，Evil Genius

今天我们分享一个分析方法思路。

在HD 16um（8um通常不采用）、Stereo-seq bin50以下（25um以下），在做单细胞空间联合分析的时候，采用的是RCTD的double模式，那么在这种模式下，细胞类型之间的空间关系该如何分析呢？

这内容本来放在2025培训班上讲的，不过既然粉丝们要求了，今天就分享出来。

其实在很早之前（2022年）进行空间培训的时候谈到过这个问题，大家可以看一下视频空间转录组高级个性化数据分析第一期

相信大家也看了很多的文献，关于double模式下的分析其实参考很少，Stereo-seq发表的文章大多也在bin60以上，采用更多的类似full模式的分析。

今天我提一篇文献，很有参考价值，很早之前分享过，不过我发现认真好好学习的人，真的很少。

其实最关键的地方，大家看一下

采用的是bin20进行分析（10-15um），接近HD的16um精度。

那么注释的方式是什么呢？简而言之，pattern模式。

就是在double的模式下，将注释分为两层，first pattern和second pattern。

那么在这种思路下，空间细胞类型的关系问题就演化成了两层pattern的细胞关系问题。

那么细胞邻域富集，细胞邻域的差异分析等均可以采用squidpy进行分析了，不过个人推荐cellcharter。

其中最为关键的是，此时的空间配受体分析直接与细胞类型相关联。

分享到这里，大家应该知道怎么分析double模式了吧。

最后写上封装好的RCTD的代码，其实脚本封装大家一定要学会。

# 创建解析器
parser <- ArgumentParser(description = "RCTD analysis")

# 添加参数
# 创建解析器
parser <- ArgumentParser(description = "RCTD analysis")

# 添加参数
parser$add_argument("-s", "--scrna", required = TRUE, help = "scrna file")
parser$add_argument("-t", "--type", required = TRUE, help = "scrna file type",choices = c('csv','rds'))
parser$add_argument("-a", "--anno",help = "scrna file type;csv")
parser$add_argument("-p", "--spatial", required = TRUE, help = "spatial rds file")
parser$add_argument("-o", "--output", required = TRUE, help = "Output directory for saving results")
parser$add_argument("-m", "--mode", choices = c("doublet", "full","multi"), required = TRUE, help = "RCTD mode")
parser$add_argument("-n", "--name", required = TRUE, help = "sample name")

scrna = args$scrna
type = args$type
anno = args$anno
spatial = args$spatial
outdir = args$output
mode = args$mode
sample = args$name

完整脚本如下：

#! usr/bin/R
####zhaoyunfei
####20241202

library(argparse)
library(Seurat)
library(ggplot2)
library(patchwork)
library(dplyr)
library(spacexr)
library(ggplot2)

# 创建解析器
parser <- ArgumentParser(description = "RCTD and CC analysis")

# 添加参数
parser$add_argument("-s", "--scrna", required = TRUE, help = "scrna file")
parser$add_argument("-t", "--type", required = TRUE, help = "scrna file type",choices = c('csv','rds'))
parser$add_argument("-a", "--anno",help = "scrna file type;csv")
parser$add_argument("-p", "--spatial", required = TRUE, help = "spatial rds file")
parser$add_argument("-o", "--output", required = TRUE, help = "Output directory for saving results")
parser$add_argument("-m", "--mode", choices = c("doublet", "full","multi"), required = TRUE, help = "RCTD mode")
parser$add_argument("-n", "--name", required = TRUE, help = "sample name")
# 解析命令行参数
args <- parser$parse_args()

scrna = args$scrna
type = args$type
anno = args$anno
spatial = args$spatial
outdir = args$output
mode = args$mode
sample = args$name

if(! file.exists(outdir)){dir.create(outdir,recursive = TRUE)}

if (type == 'csv'){

ref = read.table(scrna,header = T,row.names = 1,check.names = F,sep = '\t')

counts <- ref

if (!is.null(anno)){

anno = read.table(anno,header = T,sep = '\t',row.names = 1)

anno = anno[colnames(counts),]

cluster <- as.factor(anno$CellType)

names(cluster) <- colnames(ref)

}else{stop("错误：必须提供 --anno 参数", call. = FALSE)}

nUMI <- colSums(ref)

names(nUMI) <- colnames(ref)

reference <- Reference(counts, cluster, nUMI)

}else{ref = readRDS(scrna)

counts <- ref@assays$RNA@counts

cluster = ref$celltype

names(cluster) <- colnames(ref)

nUMI <- colSums(counts)

names(nUMI) <- colnames(ref)

reference <- Reference(counts, cluster, nUMI)

}

spatial =  Load10X_Spatial(data.dir = spatial, bin.size = c(16))

spatial <- subset(spatial, subset = nFeature_Spatial.016um > 200)

counts <- spatial@assays$Spatial.016um@layers$counts

colnames(counts) = colnames(spatial)

rownames(counts) = rownames(spatial)

coords <- GetTissueCoordinates(spatial)

coords = coords[,c(1,2)]

colnames(coords) <- c("x", "y")

coords[is.na(colnames(coords))] <- NULL

query <- SpatialRNA(coords, counts, colSums(counts))

RCTD <- create.RCTD(query, reference, max_cores = 8)

RCTD <- run.RCTD(RCTD, doublet_mode = mode)

spatial@meta.data$first_type = RCTD@results$results_df$first_type

spatial@meta.data$second_type = RCTD@results$results_df$second_type

pdf(paste(outdir,paste(sample,'RCTD.first_type.pdf',sep = '.'),sep = '/'))

SpatialDimPlot(spatial, group.by = "first_type", label = T, repel = T, label.size = 4)

dev.off()

pdf(paste(outdir,paste(sample,'RCTD.second_type.pdf',sep = '.'),sep = '/'))

SpatialDimPlot(spatial, group.by = "second_type", label = T, repel = T, label.size = 4)

dev.off()


saveRDS(spatial,file = paste(outdir,paste(sample,'HD.RCTD.rds',sep = '.'),sep = '/'))

生活很好，有你更好。

本文标签： 方法更新