脚本更新—实现visium数据的空间距离分析2

作者，Evil Genius

这一篇我们要更新一下我们的脚本，在文章中脚本更新----实现visium数据的空间距离分析中提到了针对visium对两种细胞类型的距离度量计算，由于时间限制并没有很完善，这一篇我们来完整的补充一下。

首先大家要考虑一下数据结构的问题，Seurat的V4和V5结构是不一样的，所有脚本在运行的时候要学会抽取信息。

其中坐标信息的变化，其中V4结构在

spatial_coords <- seurat_obj@images$spatial@coordinates

####或者
spatial_coords <- seurat_obj@images[[slice]]@coordinates

数据格式如下，坐标信息包含六列，分析的时候默认采用scale后的坐标

V5结构不同，坐标在

spatial_coords <- seurat_obj@images$slice1@boundaries$centroids@coords
head(spatial_coords)

还有点与点的单位距离

seurat_obj@images$slice1@boundaries$centroids@radius
130.0574

V5结构的10X HD数据坐标也在相同的位置

可以看到V5结构的空间坐标只有两列。

那么在对不同结构分析空间距离的时候，就需要进行针对性的代码分析。

library(Seurat)
library(SeuratData)
library(tibble)
library(ggplot2)
library(patchwork)
library(dplyr)
library(Rmagic)
library(reticulate)
source("citeseq_function.R")
dataspat <- readRDS("data/spatial/HCC1T.rds")

第一步，单细胞空间联合，这个就不做了，R版本用RCTD做一下就行，拿到解卷积后的矩阵，存放在meta.data下面。

第二步，中性粒细胞与肿瘤的距离分析（大家自己的项目要选择自己的目标细胞类型）。

首先获取目标细胞，这里比如取中性粒细胞的坐标（定义细胞含量大于0.9），大家看这个代码就知道单细胞空间联合后的矩阵放在了meta.data下面，针对V4、V5要抽取对应的位置。

info <- [email protected]
 <- rownames(info)[info$Neutrophil > 
                          quantile(info$Neutrophil, probs = c(0.9))]

那么接下来获取肿瘤细胞（第二种标签）的坐标,这里的ident就是形态学注释的结果，标注这里是肿瘤区域

select2 <- rownames(info)[info$ident == "tumor"]

这个地方就开始有了不一样的内容，就像图中的那样，select1在肿瘤区域内，距离是0，相邻距离就是1，隔了一个spot，距离就是2，以此类推,spot_distance是一个自定义函数，放在最后,这里注意如果距离为0，在代码上会显示select1和select2有重复。

dist <- spot_distance(dataspat = dataspat,
                      select1 = select1,
                      select2 = select2)

到这就获得了两种标签的距离矩阵信息

注意这里获得了任意两个spot之间的距离（包括自己与自己），我们首先要抽取目标spot的矩阵。

接下来我们要判断距离信息，正常情况我们想知道每个中性粒细胞距离肿瘤区域的距离信息，就要取每个中性粒spot到肿瘤区域的最短距离，也就是最小值。

input <- do.call(rbind, lapply(colnames(dist), function(i){
  data.frame(id = i, dist = min(dist[,i]))
}))

拿到这个就需要绘图了,就是图I，这时候还是单线图，表征两种标签的距离信息。

plot <- ggplot(input, aes(y=dist,x=number))+
  geom_smooth(span = 1,  method = "loess", se = TRUE) +
  theme_bw() + 
  theme(panel.grid = element_blank(),
        axis.text = element_text(colour = "black", size = 15),
        axis.title = element_text(colour = "black", size = 15),
        panel.border = element_rect(fill = NA, color = "black", linewidth = 1))

图I的另一条先是针对目标基因的变化图，这个就跟大家前面的分析相关了，画图代码是一样的。

最后来更新source的脚本，由于坐标位置的差异，大家需要针对自己的数据格式做修改，首先第一步就是拿到坐标，整理成两列，x,y，行名是barcode，下面是个示例(V5).