空间转录组学习笔记1

基本理论/背景知识参考：wx公众号：单细胞空间交响乐、空间组学；B站：七夜听雪ii的个人空间-七夜听雪ii个人主页-哔哩哔哩视频

代码参考：Seurat官方教程：使用 Seurat 对空间数据集进行分析、可视化和集成 • Seurat

https://mp.weixin.qq.com/s/_FwWt3LoecAmV8AzN3GafQ

3月19日 读入数据并整合：

# ST多样本整合分析
library(Seurat)
library(tidyverse)
library(RColorBrewer)

# 
input_data_dir <- "../space0304/raw"
sample_list <- list.files(input_data_dir)
sample_list

# 选择样本
sample_list <- sample_list[c(2:3)]
sample_list

## 文件路径
samples_dir <- sample_list %>% file.path(input_data_dir, .)
samples_dir

## 切片ID
sample_names <- c("", "")



## 依次读取每个样本
sample_objects <- purrr::map(1:length(sample_list), function(x) {
  ## 读取数据
  one_dir <- samples_dir[x]
  sample_id <- sample_list[x]
  slice_id <- sample_names[x]
  
  sample_object <- Load10X_Spatial(
    data.dir = one_dir,
    filename = "filtered_feature_bc_matrix.h5",
    assay = "Spatial",
    slice = slice_id,
    filter.matrix = TRUE
  )
  
  ## 设置样本名称和元数据
  sample_object@project.name <- sample_id
  sample_object@meta.data$orig.ident <- slice_id
  
  ## 重命名细胞ID
  sample_object <- RenameCells(object = sample_object, add.cell.id = slice_id)
  
  return(sample_object)
})



# 各样本分别进行SCT转换
sample_objects <- lapply(
  sample_objects,
  SCTransform,
  assay = "Spatial",
  method = "poisson"
)

# 合并所有样本
ST <- merge(
  sample_objects[[1]],
  y = sample_objects[2]
)

# 设置默认assay为SCT
DefaultAssay(ST) <- "SCT"

# 合并所有样本的可变特征
VariableFeatures(ST) <- c(
  VariableFeatures(sample_objects[[1]]),
  VariableFeatures(sample_objects[[2]])
)

# 运行PCA
ST <- RunPCA(
  ST,
  assay = "SCT",
  verbose = FALSE
)

# 使用harmony进行批次矫正
library(harmony)
ST <- RunHarmony(
  ST,
  reduction = "pca",
  group.by.vars = "orig.ident",
  reduction.save = "harmony"
)

# 聚类分析
ST <- FindNeighbors(
  ST,
  reduction = "harmony",
  dims = 1:30
) %>%
  FindClusters(resolution = 0.4)

# 运行UMAP
ST <- RunUMAP(
  ST,
  dims = 1:30
)

display.brewer.all()
cols <- colorRampPalette(brewer.pal(9,'Set1'))(15)
DimPlot(ST, reduction ="umap",cols = cols,group.by = "orig.ident")
SpatialPlot(ST,group.by ='seurat_clusters',ncol = 2)& scale_fill_manual(values = cols)
SpatialFeaturePlot(ST, features = c("Pecam1","Cdh5"), alpha = c(0.1, 1))

library(scRNAtoolVis)
cellRatioPlot(object = ST,
              sample.name = "orig.ident",
              celltype.name = "seurat_clusters",
              flow.curve = 0.5,
              fill.col = cols)+RotatedAxis()