Squidpy空间转录组可视化分析:从AnnData到组织切片 空间转录组数据越来越常见但很多人做完聚类、注释、UMAP 之后就把数据交给下游富集分析了组织切片上的位置信息基本没有参与后续判断。这个系列第二篇我直接从 Scanpy 读进来的 AnnData 讲起完整演示 Squidpy 在空间数据的可视化分析里能做哪些事以及每种图上到底能读出什么、读不出什么。如果你已经能跑通 Scanpy 的常规单细胞流程但还没把基因表达或细胞类型映射回真实的组织切片坐标这篇正好帮你把最后那块拼图补上顺带减少一点我当年走了很久的弯路。1. 当你已经做完聚类空间数据里真正比单细胞多出来的部分1.1 单细胞维度看不到的信号细胞邻域效应常规 scRNA-seq 的分析逻辑是把组织打碎成单细胞悬液然后按照基因表达谱给细胞分类。这个设计天然丢掉了一个关键信息细胞和谁靠在一起。肿瘤微环境里CD8 T 细胞是不是真的贴着巨噬细胞还是只是同时存在于同一块组织中这两种情况对机制解释是完全不同的。前者可能意味着 T 细胞与巨噬细胞存在物理接触形成的免疫抑制环路后者可能只是两个独立过程在同一组织里共存。表达谱聚类永远无法直接回答贴不贴这个问题只有看到空间坐标才能回答。空间转录组的核心增量就在这里除了表达矩阵还有每个捕获单位spot 或 cell的物理坐标。有了坐标我们可以计算邻居、统计邻域组成、找空间聚集的基因、判断细胞类型之间的共定位关系。这些指标不是给图加点装饰而是真正产生了新的生物学结论的维度。1.2 Visium、Slide-seq 与成像平台对可视化逻辑的影响不同平台的坐标含义差异很大但聊可视化之前必须先搞清楚点之间的邻域尺度到底是什么。10x Visium 是规则排列的捕获探针spot 之间中心距约 100 微米每个 spot 通常捕获多个细胞坐标单位是像素Slide-seq 用的是随机排布的 barcoded beads单个 bead 可能只覆盖一个到几个细胞坐标单位通常是微米MERFISH 这类成像平台本身就是单细胞分辨率的每个点就是一个确定的细胞坐标精度最高。平台差异直接决定你在 Squidpy 里怎么设参数。Visium 这种网格化排布可以用coord_typegrid走规则邻居邻居的定义很清晰Slide-seq 和 MERFISH 是散点式分布必须用coord_typegeneric设定物理半径或 K 近邻否则会把相隔很远的点硬拉成邻居。可视化时也是如此Visium 的斑点适合画成六边形点Slide-seq 画成小圆点更不容易互相遮盖。后面第 5 章我会把邻居构建的参数展开讲清楚。1.3 Scanpy 与 Squidpy 的分工边界很多刚接触空间数据的同学会问Scanpy 不是已经有sc.pl.spatial了吗为什么还要 Squidpy我的理解是Scanpy 的核心职责是表达矩阵层面的分析PCA、聚类、UMAP、marker 基因这些事它在单细胞领域已经做得很成熟空间可视化只是它的一个辅助功能参数远没有专门工具灵活。Squidpy 则是在 AnnData 这个数据结构上专门做空间分析它既提供绘图接口也提供一套完整的空间统计方法——邻域富集、共现系数、空间自相关、形态学特征等。实际项目里我的分工习惯是预处理、聚类、注释继续用 Scanpy到了聚类之后怎么办这一步切到 Squidpy。Scanpy 负责哪个细胞群是什么Squidpy 负责这些细胞群在空间里如何组织。两者共用一个 AnnData流程衔接非常顺不需要来回转换数据格式。2. 从 Space Ranger 到 AnnData把坐标、图像和表达矩阵对齐在同一张表里2.1 标准 Visium 读取流程与 Squidpy 依赖的字段Squidpy 的所有分析都建立在 AnnData 结构上但它并不是什么数据都能直接吃。关键点在于表达矩阵放在adata.Xspot 的坐标放在adata.obsm[spatial]组织图像和缩放因子放在adata.uns[spatial]。三个部分必须对齐到同一个索引Squidpy 才能把表达量原样映射到切片图上。读 10x Visium 数据时我一般直接走 Scanpy 的读取接口它能一次性把 Space Ranger 输出解析成上述结构import scanpy as sc import squidpy as sq adata sc.read_visium( path/to/10x_visium/, library_idC1, load_imagesTrue, ) adata.var_names_make_unique()这里library_id要和 Space Ranger 输出目录里的样本名一致load_imagesTrue才会把组织切片的高清图和低清图读进来。读完以后可以快速验证一下结构print(adata.obsm[spatial].shape) # (n_spots, 2)x/y 像素坐标 print(adata.uns[spatial].keys()) # 里面是 library_id - 图像和 scalefactors print(adata.uns[spatial][C1][images].keys()) # 通常有 hires 和 lowres很多非 Visium 平台的数据没有sc.read_visium这么方便的入口。如果只有一张坐标表就手动把坐标列写进adata.obsm[spatial]import pandas as pd coord_df pd.read_csv(slide_seq_coords.csv, index_col0) adata.obsm[spatial] coord_df[[x, y]].values只要列名对齐、顺序和adata.obs_names一致Squidpy 不会在意坐标来自什么平台后续的空间邻居构建和可视化都照常进行。2.2 网格数据与通用坐标数据在 spatial_neighbors 中的差异这是最容易踩坑的地方。Visium 的 spot 排布非常规矩六个邻居环绕一个中心点所以 Squidpy 针对这种结构提供了coord_typegrid的专用逻辑。网格模式下n_rings表示向外扩散多少圈邻居n_rings1就是最内层六个邻居n_rings2会把第二圈的十二个点也算进来。如果不确定该用几圈我通常先用 1 到 2 看结果稳定性再决定是否加大。Slide-seq、Stereo-seq、MERFISH 这些平台没有规则网格必须用coord_typegeneric。这时候要么指定radius单位跟着你的坐标单位走通常是微米要么指定n_neighs做 K 近邻。我个人更推荐设radius而不是死磕 K固定 K 会导致高密度区域邻居集中在非常小的物理范围低密度区域反而拉到很远生物学上可能完全不是相邻固定半径语义就清楚得多只是不同样品密度不同半径需要试几个值再看邻居数量分布。构建的时候我会在代码里写清楚sq.gr.spatial_neighbors(adata, coord_typegeneric, radius30)跑完以后随便抽几个 spot 检查一下邻居的数量如果平均邻居数少于 3说明半径设小了如果超过 20说明半径偏大图的局部特异性会被稀释。这一步不需要追求精确只要邻居数在合理范围即可。2.3 为什么要把图像和坐标放进同一个对象再画图没有图像的情况下spatial_scatter当然也能画坐标散点图但组织边界、解剖结构完全看不见点与点之间的关系很难解读。Squidpy 的imgTrue会把 AnnData 里存的组织切片图像作为底板自动使用对应的缩放因子把坐标对齐。这个设计比手动把图像和散点叠在一起更可靠因为它把 alignment 的复杂度封装进了函数内部而不是让每个分析者自己去处理像素和坐标的换算。刚开始上手的人不需要关心具体换算过程但要意识到一点千万不要自己缩放adata.obsm[spatial]来配合图像。除非你有明确理由否则坐标和图像之间的缩放由sc.read_visium读入的 scalefactors 负责Squidpy 绘图时会自动处理。手动改坐标往往导致点飘到组织外面很难找回头。3. 预处理阶段就要做的三个空间检查QC、组织覆盖与批次伪影3.1 先画 count / gene 分布再决定过滤阈值空间数据预处理最忌讳照搬单细胞阈值。单细胞样本里一个低质量细胞可能是测序深度不够但空间转录组里一个低 quality 的 spot 还可能是组织切片厚度不均、边缘区域细胞密度低或者坏死区。与其机械地把min_genes和min_counts一划不如先把这些指标画回到组织图谱上看看分布是否有空间规律。adata sq.datasets.visium_hne_1() sq.pl.spatial_scatter( adata, color[total_counts, n_genes_by_counts], imgTrue, alpha0.8, ncols2, )我自己看这类图时关注三个点第一低表达 spot 是否集中在某一侧边缘如果是多半是组织切片边缘的物理切割区域而不是测序质量问题过滤时不必太严格第二组织内部有没有零星出现的空洞有些是气泡伪影有些是探针捕获效率低的区域如果空洞面积大后续差异分析要把这些 spot 主动剔除第三整片组织的 count 中位数是否呈现明显的一侧高、一侧低梯度分布这种梯度往往意味着切片厚度不均匀或样本降解对下游表达差异的影响比几个离群 spot 更严重。所谓空间视角的 QC核心就是不让数字分布脱离组织背景去理解。单看小提琴图你只知道有离群值不加背景你会误判离群值该不该删。3.2 组织残留与边缘点怎么剔除Visium 捕获过程中探针不一定只落在组织上切片边缘、折叠、褶皱区域也可能产生一些伪 spot。它们的特征是坐标不在组织内部但表达数据依然存在。用imgTrue看图时这些点会飘在组织区域外面或者沿着边缘排成一条线。最简单的处理方法是在后续分析里过滤掉它们而不是保留下来干扰邻居关系。具体做法因项目而异。比较直接的方法是用组织图像做阈值分割找到组织覆盖区域只保留坐标落在组织内部的 spot。这一步如果你对图像处理不熟可以用手动圈选或者干脆按坐标画图后人工确定边界。不用做得特别精细空间分析对边缘少数 spot 的容忍度比单细胞分析高但要确保没有大块的非组织点在空间邻居图里被当成真实邻居否则邻域富集很容易在切片边缘产生假信号。3.3 批次效应不能只看 UMAP要回到空间图里验证单细胞批处理效应通常靠 UMAP 观察 samples 是否混在一起但空间数据的批次效应还有一个特殊表现不同样本之间的空间坐标分布完全不同即使同一个细胞类型在不同切片上的空间夹带也可能有系统差异。如果只看 UMAP 把样本混起来了就忽略了一个问题——样本 A 的高表达区域可能集中在切片左上角样本 B 的高表达区域可能均匀分布。这种差异不是简单的技术批次还可能来自样本间组织结构的真实不同。我的习惯是在 Harmony 或 BBKNN 校正之后把批次的标签直接画回空间图谱sq.pl.spatial_scatter(adata, colorbatch, imgTrue)如果同一 batch 的 spot 在空间上连成一个边角区域说明批次信息大概率来自组织切片本身的位置效应比如取样部位不同校正时就要特别小心不要把所有位置效应都当成技术批次抹掉。如果批次标签在空间图上是随机混合的但表达谱上还是能看出 split那时候再考虑更激进的技术校对比比较安全。4. 把表达量画回组织切片Squidpy 可视化函数的核心参数与组合用法4.1 单基因可视化与连续色标最简单的用法是画单个基因。我建议在color参数里传基因名然后打开图像底板sq.pl.spatial_scatter( adata, colorLYZ, imgTrue, alpha0.8, cmapmagma, )这张图能立刻告诉你 LYZ 高表达区域是集中在一个结构块还是沿着某种边界分布还是在全组织里弥散。连续色标我默认用magma或viridis它们对色觉障碍人群更友好而且颜色深浅的感知单调性比 rainbow 好得多。空间图里最忌讳用彩虹色做连续基因不同波段之间的跳跃感会让人误判高表达区域。如果你有多个基因要看直接把列表传给colorSquidpy 会自动排成多个子图sq.pl.spatial_scatter( adata, color[LYZ, CD68, CD3D, MS4A7], imgTrue, cmapmagma, ncols2, )在正式出图前我会先用ncols2或ncols3扫一遍几十个候选基因效率比一个个画高很多。扫的时候不要开savefig直接屏幕上过一遍即可。4.2 离散标签与多基因面板基因表达是连续值但聚类标签、细胞类型注释是离散值两者的视觉逻辑完全不同。离散值一定不要用连续色标会让不同组别的颜色深浅产生高低错觉。Squidpy 对离散标签会自动分配调色板但自动分配可能让相邻细胞群的色差太小肉眼看不清边界。我一般手动指定palettesq.pl.spatial_scatter( adata, colorleiden, imgTrue, alpha0.8, palettetab20, legend_locright margin, )legend_locright margin可以在图例很长时避免遮住组织主体。如果聚类数特别多比如 20 个以上的群图例会非常臃肿这时候我倾向于只画重点细胞类型比如把其他群统一设成浅灰色只突出目标群再配上groups参数限制显示范围。如果想要看细胞类型之间是否在空间上有嵌套关系直接把多个标签并排画sq.pl.spatial_scatter( adata, color[leiden, cell_type, total_counts], imgTrue, ncols3, )这种并排图是我在做报告时用得最多的形式既能快速展示聚类结果又能对照细胞类型注释还能顺带观察 count 分布是否和某个细胞类型区域重合。三张图放在一起很多空间规律一眼就能看出来。4.3 用图像做底板img、crop_coord 与 alphaimgTrue默认显示 AnnData 里的高清组织图但有些组织切片太大全图分辨率不够或者你只关心某个局部区域。这时用crop_coord截取局部分析sq.pl.spatial_scatter( adata, colorleiden, imgTrue, crop_coord[1000, 3000, 2000, 4000], )crop_coord的顺序在不同版本里可能略有差异我习惯先画全图再目测一个大概的像素范围一次不行就调两次。局部区域放大后spot 之间的接触关系和边界形态比全图清楚得多尤其适合判断两个细胞群之间是清晰分界还是互相渗透。alpha参数控制点在背景图上的透明度。我的个人经验是连续基因图可以把 alpha 放到 0.8尽量保留图像细节离散标签图如果簇之间边界很清晰alpha 0.6 就够太透明会显得颜色发白太不透明又会盖住组织纹理。这个参数没有标准答案和你的组织切片染色深浅有关建议先跑一张默认参数再上下浮动调整。5. 空间关系不漏在图上邻域富集、共现曲线和 Morans I 的组合解读5.1 构建空间邻居图是一切空间统计的前提要让 Squidpy 做空间分析第一步不是画图而是构建空间邻居图。这个过程会生成两个矩阵spatial_connectivities和spatial_distances分别记录每个 spot 和哪些 neighbor 相连、物理距离是多少。后面所有邻域富集、共现分析、空间自相关都依赖这两个矩阵所以这一步跑错了后面全是错的。Visium 样本我常用sq.gr.spatial_neighbors(adata, coord_typegrid, n_rings2)非规则平台我常用sq.gr.spatial_neighbors(adata, coord_typegeneric, radius30)跑完以后我会习惯性地看一个 summarize 信息import numpy as np n_neigh np.array(adata.obsp[spatial_connectivities].sum(axis1)).ravel() print(np.median(n_neigh), np.percentile(n_neigh, [10, 90]))如果中位数邻居数远小于预期回头检查坐标单位或者半径设定。这一步总共花不到半分钟但能避免后面大量无效计算。5.2 邻域富集热图的判读逻辑邻居图构建好以后邻域富集分析是很多空间项目的第一张统计图sq.gr.neighborhood_enrichment(adata, cluster_keyleiden) sq.pl.neighborhood_enrichment(adata, cluster_keyleiden)输出是一张热图每个格子的值是一个 Z-score代表某个细胞群在另一个细胞群周围出现的频率是否显著高于随机预期。Z-score 大于 0 意味着两个群更经常互相靠近小于 0 意味着两个群倾向于互相排斥。读这张图有个容易犯的错误只看 Z-score 的正负不看关联 p 值。Squidpy 返回的结果里同时有pval_norm和pval_norm_fdr_bh必须用 FDR 校正后的值筛选显著互作否则大样本里到处都是名义显著的结果。我的操作是先看 FDR 小于 0.05 的格子再在这些格子里挑 Z-score 绝对值最大的几对做重点解读。如果两个细胞类型确实在组织边缘不接触但在组织内部其实是贴着的邻域富集在全局水平可能检测不到显著信号因为边缘区域的非接触把效应稀释了。这时建议把边界内外分开看或者结合下一节的共现曲线看不同距离梯度上的接触情况。5.3 共现曲线距离梯度里藏着更多信息邻域富集只回答是否倾向共定位而共现分析把距离切成多个 bin看细胞类型的共现概率怎么随距离变化。Squidpy 的调用方式很简单sq.gr.co_occurrence(adata, cluster_keyleiden) sq.pl.co_occurrence(adata, cluster_keyleiden)曲线横轴是距离纵轴是条件概率。如果两条细胞类型的曲线在很近距离处就高于各自的整体比例说明它们紧密相邻如果曲线在远处才升高说明它们虽然存在于同一组织但不直接接触。我最常看的是一号和二号细胞类型在曲线起始位置的交叉起始高度很高通常意味着 spot 自身和邻居直接接触起始高度低但随后爬升则可能是两个群有共同的组织区域边界但不互相混合。用邻域富集 共现曲线两组结果互相印证比只看一张热图可靠很多。5.4 Morans I 帮你快速筛出全局聚集的基因空间自相关是另一个高频命令。Morans I 本质上在回答一个问题某个基因的表达是高高低低地随机分布还是在空间上有明显的聚集模式如果是聚集那聚集尺度是多大genes adata.var_names.tolist() moran sq.gr.spatial_autocorr( adata, modemoran, genesgenes, n_perms1000, n_jobs4, ) moran moran.sort_values(I, ascendingFalse) moran.head(10)结果里I越大代表这个基因的表达越倾向于空间聚集同时会给出 p 值和 FDR 校正后的 p 值。实际项目里我不会对成千上万个基因逐个做空间自相关一般先用var_names子集筛掉零表达占比过高的基因再跑百来个感兴趣的目标基因这样计算时间可控。还有一个容易被忽略的点Morans I 的尺度高度依赖邻居半径。同一个基因邻居半径设 30 微米和设 100 微米结果可能完全不同。半径小反映的是局部邻域模式半径大反映的是组织区域尺度模式。所以出组学报告时一定要在方法里写清楚邻居参数否则读者无法判断你的空间聚集是什么尺度上的聚集。6. 画图阶段踩过的坑坐标镜像、配色误导与大画布卡顿6.1 坐标镜像与切片方向不一致我第一次拿一批非标准输出的坐标数据画图时点全部诡异地对不到组织图像上有的左右镜像有的上下颠倒。原因很简单不同来源的坐标原点定义不同。有的平台原点在左下角有的在左上角matplotlib 默认 y 轴方向又和图像存储方向不一致。解决办法不是直接手动翻转坐标而是让 Squidpy 用imgTrue自动对齐如果对齐后仍然反了再考虑检查原始坐标文件和图像元数据。如果确实需要手动改正方向不要改adata.obsm[spatial]本身用ax.invert_yaxis()或者在画图前复制一份坐标做翻转这样不至于污染后续分析。做完之后最关键的一步是肉眼验证找一两个形态特征明显的基因或者组织结构确认图上的高表达区域和病理图像上的结构位置确实对应。6.2 配色如何影响空间结论空间图对配色的敏感度比普通散点图高得多因为观众对组织上哪里亮、哪里暗会有下意识的生理反应。连续基因我建议固定用magma或viridis并且在多图对比时保持同一套 cmap 和 vmin/vmax否则两个基因之间根本无法比较绝对值。离散标签如果簇数多自动 palette 经常出现两个相邻簇颜色太像的尴尬情况。我踩过一次很深的坑有两个细胞群在组织里其实是相互交织的但自动配色让它们在图上看起来像同一个群差点误判成细胞类型边界意外清晰。后来我改成手动检查调色板必要时把关注组之间的色差拉大。调色板这个看似是美观问题在某些情况下会直接影响分析判断。6.3 大画布与输出文件全组织切片几万个 spot 在 notebook 里的渲染压力不小。默认 dpi 情况下保存 PDF 会生成几十 MB 的文件打开都很卡。我现在的习惯是能存 PNG 就存 PNG出报告时保持 dpi200 左右足够清晰如果需要矢量图用于论文投稿就只对有代表性的局部区域用crop_coord出矢量图全图用高分辨率 PNG。还有一个小技巧spot_size不要调得过大。点的大小会直接影响视觉上这片区域是否表达高的感受过大时低表达区域的 gap 会被填满看起来好像所有区域都有表达过小时高表达区域又显得稀疏。我一般从默认值出发边看图边调直到组织结构和表达信号都能看清为止。6.4 结果解读的常见误区把统计图当成确定性证据最后说一个方法论层面的坑。邻域富集的 Z-score、共现曲线的条件概率、Morans I 的正负值这些指标都是描述性质的工具不能替代因果推断。你发现 T 细胞和巨噬细胞显著共定位只能说两者空间关系不随机不能直接说它们之间一定有互作机制你也无法从共定位推出谁主动迁到谁旁边。做空间分析要时刻记得空间关系的价值在于缩小候选机制的范围而不是直接证明机制。我写报告时会把空间统计结果称为线索生成阶段后续还需要用配体受体分析、邻域组成差异或者功能实验去验证。这样表达既不夸大结论也能让合作者明白空间分析在整个证据链里的位置。回到实际使用我现在做一套空间数据可视化分析流程基本固定读入数据后先画 QC 指标确认组织质量然后标准化聚类再看几组候选基因和细胞类型的空间分布接着构建邻居图做邻域富集和共现最后用 Morans I 筛一批空间聚集基因进入下一步分析。这里面最花时间的其实不是跑命令而是根据每张图的结果不断调整对组织的理解。空间数据最大的优势就是它的结果能落回实体组织上抓住这个感觉后面学什么新方法都快。