Skip to content

单细胞空间转录组: 基础数据分析(基于 Seurat 的单样本流程) ​

作者: SeekGene
时长: 26 分钟
字数: 6.2k 字
更新: 2026-10-13
阅读: 0 次
空间转录组 Notebooks

1. 模块简介 ​

本模块基于 Seurat V4 框架开发,专为处理和分析 SeekSpace 空间转录组单样本测序数据 而设计。 Seurat 是一套功能强大的 R 语言工具集,能够对空间转录组数据进行系统性的预处理、标准化、降维和聚类分析。针对 SeekSpace 空间转录组单样本数据,本流程主要实现以下核心分析目标:

  • 数据构建:读取 SeekSpace 平台输出的表达矩阵,并与物理空间坐标和背景图像无缝结合,构建完整的空间转录组对象。

  • 质量控制:通过计算线粒体比例、UMI 数等关键指标对数据进行质控过滤,去除低质量细胞,确保下游分析的可靠性。

  • 特征提取与降维:识别高变基因 (HVGs),并执行 PCA 和 UMAP 降维,在低维空间中揭示数据的内部结构。

  • 细胞聚类与空间映射:基于基因表达谱对细胞进行无监督聚类,并将聚类结果、特征基因的表达情况映射回组织切片的真实物理空间中,直观展示单一样本内的空间异质性。

2. 输入文件准备 ​

SeekSpace 技术的配套的基础数据处理的软件是 SeekSpace® Tools,该软件可以从测序文库中识别细胞的表达信息,同时可以识别出每个细胞在空间上的位置。

经过 SeekSpace® Tools 软件处理之后,所得到的结果文件格式和目录结构如下:

text
data/                                             # 数据主目录
├── matrix/
│   └── WTHA5_filtered_feature_bc_matrix/         # 过滤后的特征矩阵及相关图像目录
│       ├── barcodes.tsv.gz                       # 细胞条形码
│       ├── features.tsv.gz                       # 基因特征
│       ├── matrix.mtx.gz                         # 表达量稀疏矩阵
│       ├── cell_locations.tsv.gz                 # 细胞空间坐标文件
│       ├── WTHA5_aligned_DAPI.png                # 组织切片的 DAPI 染色图片
│       └── WTHA5_aligned_HE.png                  # 组织切片的 H&E 染色图片
└── rds/
    └── WTHA5.rds                                 # (可选) 包含细胞分群和注释等结果的 Seurat 对象

空间坐标文件 (cell_locations.tsv.gz) 说明:

  • 第 1 列:barcode,顺序与 matrix.mtx.gz 中的细胞顺序一致。
  • 第 2 列 & 第 3 列:分别为该 barcode 所代表的细胞在空间芯片上的像素坐标。
  • 物理距离换算:SeekSpace 技术中一个像素点的大小约为 0.2653 微米,将像素点的坐标乘以 0.2653 即可转换计算细胞在真实物理空间上的距离。
R
# --- 环境准备 ---
library(Seurat)
library(tidyverse)
library(ggplot2)
output
Attaching SeuratObject

Warning message:
“package ‘tibble’ was built under R version 4.3.3”
Warning message:
“package ‘tidyr’ was built under R version 4.3.3”
Warning message:
“package ‘dplyr’ was built under R version 4.3.3”
Warning message:
“package ‘stringr’ was built under R version 4.3.3”
Warning message:
“package ‘lubridate’ was built under R version 4.3.3”
── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✓ dplyr 1.1.4
✓ readr 2.1.5

✓ forcats 1.0.0
✓ stringr 1.5.1

✓ ggplot2 4.0.0
✓ tibble 3.2.1

✓ lubridate 1.9.4
✓ tidyr 1.3.1

✓ purrr 1.0.2

── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✗ dplyr::filter() masks stats::filter()
✗ dplyr::lag() masks stats::lag()
ℹ Use the conflicted package () to force all conflicts to become errors

3. 创建与读取 Seurat 对象 ​

在本节中,我们提供两种将 SeekSpace 数据加载到 R 环境中并构建为 Seurat 对象的方法,用户可以根据实际情况进行选择:

  • 方法一:从原始表达矩阵创建 直接读取 SeekSpace tools 输出的 filtered_feature_bc_matrix 目录下的稀疏矩阵文件(matrix.mtx.gz、barcodes.tsv.gz、features.tsv.gz)。这里使用 Seurat 常用的 Read10X 配合 CreateSeuratObject 函数即可完成从头构建。

  • 方法二:直接读取已处理的 RDS 对象 除了输出原始矩阵,seekspacetools 还会自动生成一个封装好的 Seurat 对象 .rds 文件(通常位于 rds/ 目录下)。该对象不仅包含了完整的基因表达量矩阵,还内置了空间坐标信息。您可以直接使用 readRDS 函数将其加载,快速进入下游的细胞聚类与空间可视化分析阶段。

R
# --- 数据加载与对象创建 ---
mouse_ovary.data = Read10X('./data/matrix/WTHA5_filtered_feature_bc_matrix/')
mouse_ovary = CreateSeuratObject(counts=mouse_ovary.data,project='mouse_ovary')
mouse_ovary_rds = readRDS('./data/rds/WTHA5.rds')

4. 数据质控与过滤 ​

本节对整合后的数据进行质控指标计算和可视化,然后根据质控标准过滤低质量细胞。

4.1 质控指标计算 ​

关键质控指标:

  • nCount_RNA(总 UMI 数):每个细胞检测到的总 UMI 数,反映测序深度。
  • nFeature_RNA(检测到的基因数):每个细胞检测到的基因数量,反映细胞捕获效率。
  • percent.mt(线粒体基因比例):线粒体基因表达占比,用于识别低质量或死亡细胞。

质控过滤策略(以结果为导向,无绝对标准):

  • 线粒体基因比例:过滤 percent.mt > 20% 的细胞。线粒体基因比例过高的细胞可能已死亡或受损。
  • 总 UMI 数:过滤 nCount_RNA < 200 或 nCount_RNA > 100000 的细胞。UMI 数过少表示测序深度不足,过多可能是双细胞或技术异常。
  • 检测到的基因数:过滤 nFeature_RNA < 1 或 nFeature_RNA > 10000 的细胞。基因数过少表示捕获效率低,过多可能是双细胞。

图注说明:下图展示了质控过滤前各质控指标的分布情况,按样本分组显示。

  • 横坐标:不同样本,便于比较不同样本间的质控指标差异。
  • 纵坐标:各质控指标的数值,包括 nCount_RNA(总 UMI 数)、nFeature_RNA(检测到的基因数)和 percent.mt(线粒体基因比例)。
  • 小提琴图:每个小提琴图代表一个样本,宽度表示该数值区间的细胞密度,越宽表示该区间的细胞越多。
  • 用途:用于评估原始数据质量,识别异常值,以及确定质控过滤阈值。
R
# --- 计算线粒体基因比例 ---
mouse_ovary[["percent.mito"]] = PercentageFeatureSet(mouse_ovary, pattern = "^mt-")
# --- 质控指标小提琴图可视化 ---
VlnPlot(mouse_ovary, features = c("nCount_RNA", "nFeature_RNA","percent.mito"),group.by = "orig.ident",
        ncol = 3,pt.size = 0)
output
Warning message:
“\`aes_string()\` was deprecated in ggplot2 3.0.0.
ℹ Please use tidy evaluation idioms with \`aes()\`.
ℹ See also \`vignette("ggplot2-in-packages")\` for more information.
ℹ The deprecated feature was likely used in the Seurat package.
Please report the issue at .”

4.2 质控过滤 ​

根据质控标准过滤低质量细胞。

R
# --- 质控过滤:去除低质量细胞 ---
mouse_ovary <- subset(mouse_ovary, subset = percent.mito < 20 & 
                      nCount_RNA > 200 & nCount_RNA < 100000 & 
                      nFeature_RNA > 1 & nFeature_RNA < 10000)

# 输出过滤后的细胞数
cat("过滤后细胞数:", ncol(mouse_ovary), "\n")
# --- 质控后小提琴图可视化 ---
VlnPlot(mouse_ovary, features = c("nCount_RNA", "nFeature_RNA", "percent.mito"), group.by = "orig.ident",
        ncol = 3, pt.size = 0)
output
过滤后细胞数: 67576

5. 数据预处理与降维 ​

本节对质控后的数据进行标准化、高变基因筛选、数据缩放和主成分分析(PCA)。

5.1 数据标准化与高变基因筛选 ​

数据标准化 (NormalizeData):

  • 功能:对原始 UMI 计数数据进行标准化处理,消除细胞间测序深度差异。
  • 方法:采用对数标准化方法(LogNormalize),基准测序深度为 10,000 个转录本。
  • 公式:log(1 + UMI_count / 10000)

高变基因筛选 (FindVariableFeatures):

  • 功能:识别细胞间表达变异最显著的基因,聚焦具有生物学意义的特征基因。
  • 方法:使用方差稳定变换法(vst),选取前 2,000 个高变基因。
  • 意义:高变基因通常对应细胞类型标记基因,排除技术噪音,保留真实生物学变异。

数据缩放 (ScaleData):

  • 功能:对表达矩阵进行 Z-score 标准化,为降维分析准备标准化数据。
  • 特征:每个基因的均值为 0,标准差为 1,消除基因间表达量级差异。
  • 意义:避免高表达基因主导降维结果,确保所有基因在 PCA 中平等贡献。

5.2 主成分分析 (PCA) ​

功能概述:

  • 对标准化后的高变基因进行主成分分析,提取数据的主要变异方向。
  • PCA 将高维基因表达数据降维到低维空间,保留主要生物学信号。
  • 前几个主成分通常捕获了细胞类型间的差异。
R
# --- 数据标准化 ---
mouse_ovary <- NormalizeData(mouse_ovary, normalization.method = "LogNormalize", scale.factor = 10000)

# --- 筛选高变基因 ---
mouse_ovary <- FindVariableFeatures(mouse_ovary, selection.method = "vst", nfeatures = 2000)

# --- 数据缩放(Z-score 标准化)---
# 注意:vars.to.regress 参数可用于回归掉特定变量(如线粒体基因比例、细胞周期等)
mouse_ovary <- ScaleData(mouse_ovary, features = VariableFeatures(mouse_ovary))

# --- 主成分分析 (PCA) ---
mouse_ovary <- RunPCA(mouse_ovary, features = VariableFeatures(object = mouse_ovary), reduction.name = "pca")
output
Centering and scaling data matrix

PC_ 1
Positive: Kcnma1, Mctp1, Dpf3, Tox2, Fshr, Fam13a, Prkar2b, Serpine2, Shisa6, Kcnq5
Inhbb, Prlr, Ano4, Ank3, Slc26a7, Jakmip3, Zfp385b, Inhba, Ltbp1, Stxbp6
Cables1, Chst8, Ephx2, Pgbd5, Ptprt, E330023G01Rik, Efna5, Cnmd, Lhcgr, Mro
Negative: Col4a1, Col4a2, Igfbp7, Nid1, St6galnac3, Col1a2, Col1a1, Ets1, Nhsl2, Rcsd1
Prkg1, Eng, Ltbp4, Ebf1, Shank3, Kdr, Cdh5, Col3a1, Pdgfrb, Flt1
Bicc1, Egfl7, Sparcl1, Lamb1, Rhoj, Fbn1, Mylk, Mrvi1, Frmd4b, Rbms3
PC_ 2
Positive: Cyp11a1, Runx2, Neat1, Sfrp4, Prlr, Plin4, Cemip, Ptgfr, Gm2a, Lhcgr
Fdx1, Sgk1, Mbp, Sorbs2, Idh1, Lipg, Sorbs1, Star, Scarb1, Acsl4
Tpst2, Tns3, Kcnd2, Rab27a, Myo7a, Fam129a, Abcb1b, Hsd17b7, Gsg1l, Ephb2
Negative: Csmd1, Cables1, Rad51b, Zfp385b, Kcnq5, Ank3, Mctp1, Ano4, Serpine2, Dock4
Pde7b, Kcnt2, Fam13a, Ptprd, Slc8a1, Thbs1, Fshr, Tox2, Rapgef4, Myo1e
Cfh, Prkar2b, Itih5, Slc18a2, Tmtc1, Chst8, Rbms3, Adamts19, Mki67, Slc26a7
PC_ 3
Positive: Muc16, Wnt2b, Upk3b, Plxna4, Myrf, Pkhd1l1, Upk1b, Nrg1, Ildr2, Sox6
Igfbp5, Lhx9, Pcnx2, Ppp2r2b, Hsd17b2, Cftr, Ly6e, Unc45b, Krt19, Lrrc7
Podxl, Anpep, Megf6, Grip1, Gpm6a, Stk26, Stum, Epcam, Lgr5, Krt7
Negative: Cped1, Cacna1c, Col1a2, Lama2, Ror2, Col1a1, Slit3, Cyp11a1, Frmd5, Gm17167
Egflam, Mrvi1, Pdzrn3, Dpp6, Col5a1, Dlc1, Pdgfrb, Pdgfra, Frem1, Arhgap24
Carmn, Lhfp, Prlr, 4930578G10Rik, Adcy5, Pcsk5, B130024G19Rik, Runx2, Cntfr, Col3a1
PC_ 4
Positive: Egfl7, Flt1, Shank3, Prkch, Fli1, Nos3, Pecam1, Cd93, Rapgef5, Tie1
Afap1l1, Mecom, Stab1, Cyyr1, Cdh5, Arap3, AC137513.1, Eng, Arhgef15, Pcdh17
Tek, Ptprb, Slco2a1, Plxnd1, Ushbp1, Notch4, Emcn, Nova2, Zfp366, Rasip1
Negative: Col6a2, Adgrl3, Bicc1, Muc16, Nrg1, Ildr2, Upk3b, Wnt2b, Lhx9, Myrf
Col6a1, Cacnb2, Upk1b, Pcnx2, Pkhd1l1, Sox6, Kcnab1, Cdh11, Lgr5, Ppp2r2b
Hsd17b2, Cald1, Palld, Col3a1, Krt19, Ltbp4, Cftr, Sorcs2, Ltbp2, Pdzrn3
PC_ 5
Positive: Ltbp1, Ephx2, Stxbp6, Sfrp4, Ptgfr, Shisa6, Plin4, Dpysl4, Cemip, Shroom3
Heg1, Jakmip3, Scarb1, Lipg, Lrp8, Dpf3, Sgk1, Rnf152, Hsd17b7, Gm2a
Sorbs1, Tns3, Ano4, Tubb5, Spns2, Runx2, Galnt18, Lrrn1, Vcan, Tox2
Negative: Dock2, Tbxas1, Myo1f, Ptprc, Fyb, Slc9a9, Cyth4, Cd84, Csf1r, Apobec1
Aoah, Arhgap25, Arhgap15, Dock10, Inpp5d, Wdfy4, Blnk, Lilrb4a, Ikzf1, Tmcc3
Hpgds, Ctss, Sirpa, Cd74, Pik3r5, Mrc1, Siglec1, Adap2, Prkcb, Ms4a7

6. 降维与细胞聚类 ​

本节使用 PCA 降维结果执行 UMAP 降维计算、构建邻域图并进行细胞聚类(无监督分群)。这些计算步骤是后续所有可视化展示的基础。

6.1 UMAP 降维计算 ​

功能概述:

  • 将高维的基因表达数据投影到 2D 空间,保留数据拓扑结构,为直观展示细胞间的关系做准备。

参数说明:

  • dims = 1:30:使用前 30 个主成分 (PCA) 作为输入特征。
  • reduction.name = "umap":将降维结果命名为 "umap" 并保存在 Seurat 对象中。

6.2 构建邻域图与细胞聚类 ​

构建邻域图:

  • 功能:基于 PCA 降维空间计算细胞间的相似性(如欧氏距离),构建 K 近邻 (KNN) 图,用于后续的图聚类分析。
  • 参数:使用前 30 个主成分(dims = 1:30)来定义细胞间的距离。

细胞聚类:

  • 功能:应用 Louvain 或 Leiden 算法对细胞邻域图进行社区发现(Community Detection),将细胞划分为不同的簇 (Clusters)。
  • 参数:resolution = 1,控制聚类的细致程度。分辨率值越大,划分的簇越多、越精细;反之则簇更少、更宽泛。可根据实际细胞类型的复杂度进行调整。
R
# --- UMAP 降维计算 ---
mouse_ovary <- RunUMAP(mouse_ovary, dims = 1:30, reduction.name = "umap")

# --- 构建邻域图 ---
mouse_ovary <- FindNeighbors(mouse_ovary, dims = 1:30)

# --- 细胞聚类 ---
mouse_ovary <- FindClusters(mouse_ovary, resolution = 1)
output
Warning message:
“The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
This message will be shown once per session”
10:45:57 UMAP embedding parameters a = 0.9922 b = 1.112

10:45:57 Read 67576 rows and found 30 numeric columns

10:45:57 Using Annoy for neighbor search, n_neighbors = 30

10:45:57 Building Annoy index with metric = cosine, n_trees = 50

0% 10 20 30 40 50 60 70 80 90 100%

[----|----|----|----|----|----|----|----|----|----|

*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
*
|

10:46:03 Writing NN index file to temp file /tmp/RtmphxgDME/file36273fdc955

10:46:03 Searching Annoy index using 1 thread, search_k = 3000

10:46:25 Annoy recall = 100%

10:46:25 Commencing smooth kNN distance calibration using 1 thread
with target n_neighbors = 30

10:46:26 Initializing from normalized Laplacian + noise (using irlba)

10:46:32 Commencing optimization for 200 epochs, with 3081284 positive edges

10:47:02 Optimization finished

Computing nearest neighbor graph

Computing SNN



Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck

Number of nodes: 67576
Number of edges: 2338906

Running Louvain algorithm...n Maximum modularity in 10 random starts: 0.8918
Number of communities: 28
Elapsed time: 22 seconds

7. 结果可视化 ​

数据整合分析完成后,我们需要通过多种可视化方法来展示分析结果。本节将介绍三种主要的可视化方式:1. UMAP 降维可视化:展示细胞在降维空间中的分布,便于观察不同样本和细胞群体的关系。2. 空间位置可视化:展示细胞在组织样本中的实际空间位置,用于分析空间分布模式。3. 自定义空间绘图函数:提供更灵活的可视化方式,可以结合组织图像背景进行展示。

7.1 UMAP 降维可视化 ​

功能概述:UMAP(Uniform Manifold Approximation and Projection)降维可视化是单细胞和空间转录组分析中最常用的可视化方法之一。它能够将高维的基因表达数据投影到二维空间,保留细胞间的相似性关系。

可视化内容:

  • 按聚类分组:展示不同细胞群体的分布,用于识别细胞类型和亚群。
  • 按细胞类型分组:通过组织细胞类型特异性 marker 进行区分。

图注说明(UMAP 降维可视化图:下图展示了细胞在降维空间中的分布情况,按样本或聚类分组显示。

  • 点的含义:每个点代表一个细胞。
  • 距离关系:距离相近的点表示细胞表达谱相似。
  • 颜色分组:相同颜色或标记的点属于同一类别。
R
# --- UMAP 降维计算 ---
# options(): 设置图形输出的尺寸(高度和宽度,单位为英寸)
# repr.plot.height 和 repr.plot.width: 控制 Jupyter notebook 中图形显示的大小
options(repr.plot.height=8, repr.plot.width=12)

MYCOLOR=c("#6394ce", "#2a4c87", "#eed500", "#ed5858",
          "#f6cbc2", "#f5a2a2", "#3ca676", "#6cc9d8",
          "#ef4db0", "#992269", "#bcb34a", "#74acf3",
          "#3e275b", "#fbec7e", "#ec4d3d", "#ee807e",
          "#f7bdb5", "#dbdde6", "#f591e1", "#51678c",
          "#2fbcd3", "#80cfc3", "#fbefd1", "#edb8b5",
          "#5678a8", "#2fb290", "#a6b5cd", "#90d1c1",
          "#a4e0ea", "#837fd3", "#5dce8b", "#c5cdd9",
          "#f9e2d6", "#c64ea4", "#b2dfd6", "#dbdfe7",
          "#dff2ec", "#cce8f3", "#e74d51", "#f7c9c4",
          "#f29c81", "#c9e6e0", "#c1c5de", "#750000"
          )

DimPlot(mouse_ovary, pt.size = 1,cols = MYCOLOR,group.by=c("seurat_clusters"))+
labs(x = "UMAP1", y = "UMAP2",title = "Cluster") +
  theme(panel.border = element_rect(fill=NA,color="black", size=1, linetype="solid"),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank(),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank())

7.2 细胞类型注释 ​

无监督聚类仅能将细胞分为不同的簇(Clusters,如 0, 1, 2 等),为了赋予这些簇真实的生物学意义,我们需要结合先验生物学知识(如特定细胞类型的特征 Marker 基因表达情况)为每个簇进行命名,即细胞类型注释。

R
# --- 细胞类型注释 (Cluster 映射) ---
cluster2celltype <- c("0"="theca cells",
                      "1"="theca cells", 
                      "2"="Granulosa cells", 
                      "3"="Proliferating cells", 
                      "4"="Granulosa cells", 
                      "5"="Stromal cells",
                      "6"="luteal cells", 
                      "7"="Granulosa cells", 
                      "8"="Granulosa cells",
                      "9"="Granulosa cells",
                      "10"="theca cells",
                      "11"="luteal cells",
                      "12"="Granulosa cells",
                      "13"="Endothelial cells",
                      "14"="Stromal cells",
                      "15"="luteal cells",
                      "16"="Epithelial cells",
                      "17"="Stromal cells",
                      "18"="Immune cells",
                      "19"="Endothelial cells",
                      "20"="Stromal cells",
                      "21"="Proliferating cells",
                      "22"="luteal cells",
                      "23"="Granulosa cells",
                      "24"="Epithelial cells",
                      "25"="Immune cells",
                      "26"="luteal cells",
                      "27"="Oocytes")
mouse_ovary[['cell_type']] = unname(cluster2celltype[mouse_ovary@meta.data$seurat_clusters])
R
options(repr.plot.height=8, repr.plot.width=12)
MYCOLOR=c("#6394ce", "#2a4c87", "#eed500", "#ed5858",
          "#f6cbc2", "#f5a2a2", "#3ca676", "#6cc9d8",
          "#ef4db0", "#992269", "#bcb34a", "#74acf3",
          "#3e275b", "#fbec7e", "#ec4d3d", "#ee807e",
          "#f7bdb5", "#dbdde6", "#f591e1", "#51678c",
          "#2fbcd3", "#80cfc3", "#fbefd1", "#edb8b5",
          "#5678a8", "#2fb290", "#a6b5cd", "#90d1c1",
          "#a4e0ea", "#837fd3", "#5dce8b", "#c5cdd9",
          "#f9e2d6", "#c64ea4", "#b2dfd6", "#dbdfe7",
          "#dff2ec", "#cce8f3", "#e74d51", "#f7c9c4",
          "#f29c81", "#c9e6e0", "#c1c5de", "#750000"
          )

DimPlot(mouse_ovary, pt.size = 1,cols = MYCOLOR,group.by=c("cell_type"))+
#NoLegend()+
labs(x = "UMAP1", y = "UMAP2",title = "Celltype") +
  theme(panel.border = element_rect(fill=NA,color="black", size=1, linetype="solid"),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank(),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank())

8. 添加空间坐标 ​

为了在真实的物理空间中展示细胞的分布特征,我们需要将 seekspacetools 输出的细胞空间坐标(cell_locations.tsv.gz)整合到 Seurat 对象中。

技术实现与注意事项:

  • 坐标存储方式:在此流程中,我们将二维空间坐标(X, Y)作为一种特殊的“降维结果”,利用 CreateDimReducObject 函数将其存储在 Seurat 对象的 @reductions$spatial 槽中。这样,在后续的可视化环节就可以直接利用 DimPlot 等内置函数绘制细胞的空间分布图。
  • 细胞顺序对齐(关键):构建空间坐标对象时,读入的空间坐标矩阵中细胞条形码(Barcode)的顺序,必须与当前 Seurat 对象元数据(@meta.data)中的细胞顺序严格一致。若顺序发生错乱,会导致基因表达谱与物理位置发生错误映射。
R
# --- 添加空间坐标信息 ---
spatial_df <- read.table('../SGS25030607_mouse_ovary_demo/data/matrix/WTHA5_filtered_feature_bc_matrix/cell_locations.tsv.gz', row.names = 1, sep = '\t',header = T)
colnames(spatial_df) <- c("spatial_1","spatial_2")
spatial_matrix <- as.matrix(spatial_df)
spatial_matrix_sorted <- spatial_matrix[match(row.names(mouse_ovary@meta.data),row.names(spatial_matrix)), ]
mouse_ovary@reductions$spatial <- CreateDimReducObject(embeddings = spatial_matrix_sorted, key='spatial_', assay='RNA')

9. 添加组织图像信息 ​

为了更直观地将细胞聚类结果和基因表达特征与真实的组织学形态结合起来,我们还需要将切片的病理染色图像(如 H&E 或 DAPI 荧光图像)整合到当前的 Seurat 对象中。

我们将图像本身及其物理尺寸参数(X 轴和 Y 轴的最大像素值)统一存储在 Seurat 对象的 @misc$info 槽中。这样在后续的空间绘图中,系统能够直接调用这些信息,确保细胞的坐标点精准地落在背景图像对应的解剖学结构上。

R
# --- 图像参数配置 ---
samplename = 'A5'
size_x = 55050
size_y = 19906
mouse_ovary@misc$info[[`samplename`]]$size_x = as.integer(size_x)
mouse_ovary@misc$info[[`samplename`]]$size_y = as.integer(size_y)

# --- 载入并附加 DAPI/HE 图像 ---
img = '../SGS25030607_mouse_ovary_demo/data/matrix/WTHA5_filtered_feature_bc_matrix/WTHA5_aligned_DAPI.png'  

#base64格式
img_64 = base64enc::dataURI(file = img)
mouse_ovary@misc$info[[`samplename`]]$img = img_64
#png格式
img_gg <- png::readPNG(img)
img_grob <- grid::rasterGrob(img_gg, interpolate = FALSE, width = grid::unit(1,"npc"), height = grid::unit(1, "npc"))
mouse_ovary@misc$info[[`samplename`]]$img_gg = img_grob

10. 空间绘图展示 ​

至此,Seurat 对象已经整合了基因表达谱、细胞聚类结果(Clusters)、二维降维坐标(UMAP)以及真实的物理空间坐标(Spatial)。

接下来,我们可以直接利用 Seurat 丰富的内置绘图函数(如 DimPlot、FeaturePlot 等),将分子层面的异质性特征还原到组织切片的真实解剖学空间中,从而直观地探索细胞类型的空间分布和基因的空间表达模式。

10.1 空间坐标绘图 ​

除了从原始的表达矩阵(Matrix)开始构建对象并进行全流程分析外,如果您的 data/rds/ 目录下已经存在处理好且保存了空间坐标的 Seurat .rds 文件(例如包含细胞分群和注释信息的对象),您也可以直接读取该 RDS 文件来进行空间分布的可视化分析,从而跳过前面的预处理和降维聚类步骤。

R
# --- 空间坐标水平上的聚类图 ---
options(repr.plot.height=7, repr.plot.width=15)
DimPlot(mouse_ovary, reduction = 'spatial',pt.size = 1)
R
# --- 空间坐标水平上的特征表达图 ---
options(repr.plot.height=5, repr.plot.width=20)
FeaturePlot(mouse_ovary, reduction = 'spatial', features = c('Cdh5','Nckap5'),pt.size = 0.7)

下面的代码示例将展示如何读取已有的 RDS 对象(以 WTHA5.rds 为例)并直接进行空间绘图展示:

R
# --- 空间坐标水平上的聚类图 (基于读取的 RDS 对象) ---
options(repr.plot.height=7, repr.plot.width=15)
DimPlot(mouse_ovary_rds, reduction = 'spatial',pt.size = 1,group.by = "seurat_clusters")
R
# --- 空间坐标水平上的特征表达图 (基于读取的 RDS 对象) ---
options(repr.plot.height=5, repr.plot.width=20)
FeaturePlot(mouse_ovary_rds, reduction = 'spatial', features = c('Cdh5','Nckap5'),pt.size = 0.7)

10.2 带有背景图像的空间特征图 ​

自定义了 ImageSpacePlot 和 FeatureSpacePlot 两个函数绘制细胞在空间上的分群信息和细胞连续变量的指标

R
######################绘图细胞分群函数
ImageSpacePlot = function(obj, group_by, type="DAPI", sample=names(obj@misc$info)[1], size=1, alpha=1,color=MYCOLOR){
    MYCOLOR=c("#6394ce", "#2a4c87", "#eed500", "#ed5858",
              "#f6cbc2", "#f5a2a2", "#3ca676", "#6cc9d8",
              "#ef4db0", "#992269", "#bcb34a", "#74acf3",
              "#3e275b", "#fbec7e", "#ec4d3d", "#ee807e",
              "#f7bdb5", "#dbdde6", "#f591e1", "#51678c",
              "#2fbcd3", "#80cfc3", "#fbefd1", "#edb8b5",
              "#5678a8", "#2fb290", "#a6b5cd", "#90d1c1",
              "#a4e0ea", "#837fd3", "#5dce8b", "#c5cdd9",
              "#f9e2d6", "#c64ea4", "#b2dfd6", "#dbdfe7",
              "#dff2ec", "#cce8f3", "#e74d51", "#f7c9c4",
              "#f29c81", "#c9e6e0", "#c1c5de", "#750000")
        
        raster_type <- switch(type,
                          HE = "img_he_gg",
                          DAPI = "img_gg",
                          stop("Invalid type. Must be 'HE' or 'DAPI'.")
                             )

        spatial_coord1 <- as.data.frame(obj[[group_by]])
        colnames(spatial_coord1) <- group_by
        spatial_coord2 <- as.data.frame(obj@reductions$spatial@cell.embeddings)
        spatial_coord <-cbind(spatial_coord2,spatial_coord1)
 
    
        ImageSpacePlot <- ggplot2::ggplot() + ggplot2::annotation_custom(grob = obj@misc$info[[sample]][[raster_type]],
        xmin = 0, xmax = obj@misc$info[[sample]]$size_x, 
        ymin = 0, ymax = obj@misc$info[[sample]]$size_y) +
        ggplot2::geom_point(data = spatial_coord, ggplot2::aes(x = spatial_1,y = spatial_2, color = !!sym(group_by), 
                            fill = !!sym(group_by)), size=size, alpha=alpha)+
        labs(size = group_by) + guides(alpha = "none")+ 
        ggplot2::theme_classic()+
        scale_color_manual(values = color)+ coord_fixed()
    return(ImageSpacePlot)
}

###################绘图基因表达函数
FeatureSpacePlot = function(obj, feature, type="DAPI", sample=names(obj@misc$info)[1], size=1, alpha=c(1,1),color=c("lightgrey","blue")){
    raster_type <- switch(type,
                          HE = "img_he_gg",
                          DAPI = "img_gg",
                          stop("Invalid type. Must be 'HE' or 'DAPI'.")
                             )

    spatial_coord1 <- as.data.frame(obj@reductions$spatial@cell.embeddings)
    spatial_coord2 <- FetchData(obj,feature)
    colnames(spatial_coord2) <- feature
    spatial_coord <-cbind(spatial_coord1,spatial_coord2)

    FeatureSpacePlot <-ggplot2::ggplot() + ggplot2::annotation_custom(grob = obj@misc$info[[sample]][[raster_type]],
        xmin = 0, xmax = obj@misc$info[[sample]]$size_x, ymin = 0, ymax = obj@misc$info[[sample]]$size_y) +
        ggplot2::geom_point(data = spatial_coord, ggplot2::aes(x = spatial_1, y = spatial_2,color = !!sym(feature),alpha = !!sym(feature)), size=size)+
        labs(color = feature)+
        guides(alpha = "none")+
        ggplot2::theme_classic()+
        ggplot2::scale_alpha_continuous(range=alpha)+
        scale_color_gradient(low=color[1],high = color[2])+ coord_fixed()
    return(FeatureSpacePlot)
    }
R
# --- 带有 DAPI/HE 背景的细胞分群图 ---
# DAPI/HE背景的细胞分群图
options(repr.plot.height=7, repr.plot.width=15)
ImageSpacePlot(obj=mouse_ovary, group_by = "cell_type",type="DAPI",size=0.7)
output
Ignoring unknown labels:
• size : "cell_type"
R
# --- 带有 DAPI/HE 背景的基因表达图 ---
# DAPI/HE背景的基因表达图
options(repr.plot.height=7, repr.plot.width=15)
FeatureSpacePlot(obj=mouse_ovary, feature="Cdh5",type="DAPI")

下面的代码示例将展示如何读取已有的 RDS 对象(以 WTHA5.rds 为例)并直接进行空间绘图展示:

R
# --- 带有 DAPI/HE 背景的细胞分群图 (基于读取的 RDS 对象) ---
# DAPI/HE背景的细胞分群图
options(repr.plot.height=7, repr.plot.width=15)
ImageSpacePlot(obj=mouse_ovary_rds, group_by = "seurat_clusters",type="DAPI",size=0.7)
output
Ignoring unknown labels:
• size : "seurat_clusters"
R
# --- 带有 DAPI/HE 背景的基因表达图 (基于读取的 RDS 对象) ---
# DAPI/HE背景的基因表达图
options(repr.plot.height=7, repr.plot.width=15)
FeatureSpacePlot(obj=mouse_ovary_rds, feature="Cdh5",type="DAPI")
R
0 条评论·0 条回复