十张10x Visium切片的空间转录组Python预处理与可视化实战-北京尧图网络

拿到十张 10x Visium 空间转录组切片的数据意味着什么意味着你手里的不是简单的表达矩阵而是十份“长在组织原位上的基因表达地图”。预处理做到位后面聚类、细胞类型注释、空间共定位分析全都顺预处理糊弄过去后面每一步都在还债。这篇博文我就以这十张 Visium ST 切片为例完整复盘一套基于 Python 的预处理与可视化流程从原始输出文件到最终可以直接拿去发文章的空间分布图每个环节的参数选择逻辑都展开讲。1. 先搞清楚 Visium ST 数据的“原生面貌”飞速理解十张切片意味着什么很多人拿到数据就开始跑代码但我建议你先花半天时间把数据来源和格式彻底搞清楚。空间转录组和普通单细胞测序最大的区别就是每个捕获点spot不仅有一串基因表达量还挂着它在组织切片上的物理坐标。坐标系统一旦理解错了后面画图全是乱的。1.1 从组织切片到测序数据数据是怎么一步步变出来的10x Visium 的实验原理简单说就是把新鲜冷冻或石蜡包埋的组织切片放在带有捕获探针的芯片上组织透化后 mRNA 释放被芯片上特定位置的探针捕获每个位置对应一个带有空间条形码spatial barcode的捕获点。逆转录、建库、测序之后下机的数据通过官方流程比对定量最终产出两个维度的核心信息哪个 spot 上捕获到了哪些基因、各有多少转录本表达矩阵这些 spot 在组织切片图像上的物理位置坐标文件和缩放因子 scalefactors。我在第一次处理 Visium 数据时最直观的感受是这比普通的 10x 单细胞数据多了一个“图层”的概念。你不仅在看细胞群体还在看这些群体在组织里待的“位置”。1.2 原始输出文件的“解剖图”哪些文件是后续分析的命根子数据量上来之后十张切片对应的文件结构很容易让人头大。每一张 Visium 切片经过官方流程Space Ranger处理之后核心输出包括文件/目录内容预处理阶段的作用filtered_feature_bc_matrix/过滤后的表达矩阵直接作为 Python 分析的输入spatial/tissue_positions_list.csv每个 spot 的像素坐标和行列号空间可视化时图层对齐spatial/scalefactors_json.json图像缩放参数坐标与组织图像配准spatial/tissue_hires_image.png高分辨率组织染色图可视化的底色背景spatial/tissue_lowres_image.png低分辨率组织染色图快速预览用说实话在一开始我对tissue_positions_list.csv里的坐标字段印象不深导致第一次画图时 spot 位置和图像上的组织形态完全错位那种挫败感记忆犹新。后来才意识到这个文件里的px_col、px_row是像素坐标而聚类图里用的坐标体系虽然也是基于像素但必须经过scalefactors的换算或者直接用AnnData里现成的spatial字段才能保证点是准确落在组织图像上的。1.3 十张切片的数据规模画像预处理前先摸清变量做批处理前一定要把数据的变量摸清。比如我这十张切片来自同一个个体的不同组织区域但有的切片组织形态差异很大、有的切片坏死区域多、有的透化程度不均匀。直接合并跑同一个预处理参数必然出问题。我通常拿到每张切片先做一张“质控画像”表记录总 UMIs、中位基因数、spot 数量以及组织覆盖率后续所有过滤阈值的确定都基于这份画像。2. 工具链选型为什么我把预处理和可视化全压在 Python 生态里Visium 数据预处理这条路说实话有不止一种走法。官方 Space Ranger 的 web summary 能出一堆 QC 报告但它是“观察窗口”不是“手术台”。我的选择是从 Space Ranger 的输出接管数据用 Python 生态scanpyAnnData完成从过滤到可视化的全流程。2.1 Python 生态处理空间组学数据的核心优势这里的核心优势不在于某个库有多强而在于分析路径的连续性。普通单细胞分析里scanpy 的AnnData对象已经是事实标准空间转录组数据无非是在这个标准对象里增加了obsm[spatial]和uns[spatial]两个字段。这意味着你可以在同一套代码里完成 QC、归一化、聚类、空间可视化不需要在 R 和 Python 之间来回倒数据。再叠加matplotlib强大的底层绘图控制能力几乎任何你想呈现的空间分布效果都能实现。2.2 这套流程的定位在 Space Ranger 输出之后接管所有数据操作Space Ranger 输出的filtered_feature_bc_matrix里已经做过一轮基本的 spot 过滤比如组织覆盖区域外的 spot 会被去除但离“干净可用”还有距离。这套流程的意义就是把这中间的地带填上空 spot 剔除、低质量 spot 过滤、线粒体比例过滤、深度归一化、高变基因筛选、聚类最后到空间可视化一气呵成。提示如果你打算处理超高深度数据或超大切片面板也可以考虑 GPU 加速或分布式方案但绝大多数实验室的数据规模单机 Python 完全够用。没必要一上来就上重武器。2.3 环境配置与依赖版本这一步踩的坑最少但最烦直接给一份我实际验证过的环境清单如果你用的版本和这份差距很大建议单独建环境conda create -n visum_st python3.10 conda activate visum_st pip install scanpy1.10.1 pip install anndata0.10.7 pip install squidpy1.4.1 pip install matplotlib3.8.4 pip install pandas2.2.2 numpy1.26.4 pip install scikit-learn1.4.2 pip install leidenalg0.10.2为什么单列squidpy因为它对空间邻域分析、空间统计有现成实现后面做空间自相关基因筛选时非常管用。如果你不想装太多包scanpy本身也够画基础空间图但squidpy的相关函数能省掉大量手写代码我还是建议装上。3. 预处理流程逐段拆解从 QC 阈值到聚类输入每个参数都不白给预处理不是“跑一遍流程就完事”。每一个过滤操作的背后都对应着一个生物学判断你相信哪些 spot 是真实的信号哪些 spot 是噪声。这一节我把十张切片处理过程中的每一步包括参数怎么定的、改参数时踩了什么坑都摊开讲。3.1 质量指标计算counts、genes、线粒体比例一个都不能少对每个 spot我计算三个核心 QC 指标total_counts该 spot 捕获到的总 UMI 数衡量测序深度n_genes_by_counts该 spot 检测到的基因数衡量捕获多样性pct_counts_mt线粒体基因 UMI 占比衡量细胞状态和透化质量。代码很简单import scanpy as sc import pandas as pd adata sc.read_10x_h5(filtered_feature_bc_matrix.h5) # 计算质控指标 adata.var[mt] adata.var_names.str.startswith(MT-) sc.pp.calculate_qc_metrics(adata, qc_vars[mt], inplaceTrue) print(adata.obs[[total_counts, n_genes_by_counts, pct_counts_mt]].describe())有了这张描述性统计表之后不要急着定阈值先画分布图。我通常把十张切片全部读成一个AnnData列表每张单独看分布而不是直接合并跑sc.pp.filter_cells。3.2 过滤阈值怎么定不能照搬单细胞项目的“老经验”单细胞分析里常见的n_genes 200、pct_counts_mt 20这套经验值放在 Visium 数据上不一定对。为什么因为每个 spot 捕获到的 RNA 不一定是单一细胞的转录组可能是 1 到几个细胞的总和。尤其是透化时间长导致 RNA 扩散时spot 内信号可能更“稀薄”但更混杂。我这次项目的做法是先看每张切片的分布直方图确定合理的下界和上界。以其中一张切片为例它的total_counts下四分位数是 1200中位数是 3800但有个组织覆盖不全的区域 spot 只有 150 左右明显属于边缘噪声直接剔除。线粒体比例这快也要分开看有的组织坏死区域线粒体比例能飙到 60%这种 spot 留着会严重干扰后续聚类不能手软。# 示例逐切片过滤 for name, ad in adata_list: sc.pp.filter_cells(ad, min_counts500) sc.pp.filter_cells(ad, min_genes200) ad ad[ad.obs[pct_counts_mt] 30, :].copy() adata_list[name] ad3.3 归一化与高变基因选择空间数据与单细胞分析的细微差别过滤完下一步是归一化。这一步我用的是sc.pp.normalize_total加sc.pp.log1p把每个 spot 的计数归一化到中位数或者固定值。空间转录组和单细胞在这里的区别是不同位置的 spot 捕获效率差异可能不是来自细胞类型而是来自组织透化不均匀。所以归一化目标值我强烈建议看一下最深的切片能不能和其他切片对齐否则合并分析时会有很强的切片间批次效应。高变基因选择上sc.pp.highly_variable_genes的默认参数适用于单细胞但空间转录组的 spot 数通常比单细胞少我一般把n_top_genes设为 2000-3000并且每一张切片单独算高变基因取交集后再用于后续 PCA。这样能减少切片间技术差异对高变基因列表的干扰。3.4 批次效应与多切片合并的处理思路十张切片怎么合并是预处理阶段最需要花心思的设计。如果切片来自不同个体的同一类组织合并前最好做 Harmony 或sc.pp.combat批次校正如果切片来自同一个个体的不同切面我倾向于先合并再用 Harmony 做轻量校正。这一步我会在 PCA 之后、邻域图构建之前处理因为批次效应如果不校正空间聚类很可能会先按“切片编号”分群而不是按“生物学状态”分群。import scanpy as sc adata_all sc.concat(adata_list, joinouter, labelsample, index_unique-) sc.pp.normalize_total(adata_all, target_sum1e4) sc.pp.log1p(adata_all) sc.pp.highly_variable_genes(adata_all, n_top_genes2000) sc.tl.pca(adata_all, n_comps30) # 批次校正可选 import harmonypy sc.external.pp.harmony_integrate(adata_all, keysample)注意Harmony 的key必须是adata.obs里的切片标签列名。如果你不想引入额外依赖sc.pp.combat也可以但速度慢很多。4. 空间可视化实战把表达矩阵“放回”组织切片坐标预处理做得再漂亮最终面对审稿人和合作者时还是要靠图说话。空间转录组的可视化本质上就一句话把表达矩阵的数值映射到组织切片图像的像素坐标上。这一节是实操味最重的部分从坐标系统到图层叠加逐层展开。4.1 坐标系与图像对齐理解 scalefactors 和 spot 位置的底层逻辑Visium 的数据坐标体系有两种说法一种是原始像素坐标px_col、px_row一种是组织图像缩放后的坐标scaled。在 Python 中AnnData的空间信息约定放在adata.obsm[spatial]同时adata.uns[spatial]里记录着组织图像路径和scalefactors。我之前犯过一个错误直接把tissue_positions_list.csv里的原始像素坐标塞进obsm[spatial]结果画出的聚类图所有点都缩在组织图像的左上角一小块区域。问题出在哪出在原始像素坐标是全分辨率图像的坐标而 matplotlib 显示时用的是低分辨率图像两者差了scalefactors_json.json里的tissue_lowres_scalef倍。如果你用scanpy.pl.spatial画图它会自动处理坐标系但如果你用 matplotlib 自己画就一定要先做缩放import json import numpy as np import matplotlib.pyplot as plt from PIL import Image scale json.load(open(spatial/scalefactors_json.json))[tissue_lowres_scalef] coords adata.obsm[spatial] * scale # 换算到低分辨率图像坐标4.2 基础可视化总表达量、聚类结果的空间映射先用最基础的图建立“组织形态投影感”。我用三张图做第一轮检查总表达量total_counts的空间分布检测基因数n_genes的空间分布聚类的空间分布。第一轮就能看出很多门道比如组织边缘是否有异常高表达区域边缘坏死或者透化残留、聚类结果是否在解剖结构上有清晰边界等。sc.pl.spatial(adata, colortotal_counts, spot_size1.8, showFalse) sc.pl.spatial(adata, colorleiden, spot_size1.8, showFalse, palettetab20)spot_size这个参数很多人忽略但它在空间图中非常关键。默认值在低分辨率图像上可能让 spot 覆盖面积过大或过小导致看不出组织边界。我的经验低分辨率图像一般spot_size设 1.5-2.5如果你的切片特别大、spot 数量特别多适当调小否则颜色糊成一片。4.3 单基因/自定义基因集表达的空间展示文章里最常出现的一类图就是目标基因的空间表达分布。实现上只需要把color换成基因名sc.pl.spatial(adata, colorTP53, spot_size2.0, cmapmagma, showFalse)基因表达图有两个需要留意的细节。第一如果你用默认的viridis或plasma色系低表达区域的 spot 会被淡色淹没很难和高表达区域区分。我更喜欢magma或inferno高表达亮黄低表达深紫对比更强烈。第二默认的 colorbar 是连续映射但你也可以配合阈值把低表达归零突出“有信号”的区域expr adata[:, TP53].X.toarray().flatten() expr[expr 1] 0 # 把低于阈值的表达视为无信号 adata.obs[TP53_detected] expr sc.pl.spatial(adata, colorTP53_detected, spot_size2.0, cmapmagma)4.4 进阶可视化细胞类型概率叠加与多切片对比聚类不是终点真正的空间可视化价值在于“把生物学状态翻译成位置信息”。在我们的项目里后续对接了细胞类型注释给每个 spot 打上了注释标签然后画了细胞类型分布图。更有意思的是把每种细胞类型的概率分数例如通过某个基因集打分映射到空间上形成“概率热图”能非常直观地看出不同细胞类型在组织里的分层结构。对于十张切片的对比我喜欢用matplotlib的subplot把同一基因在不同切片上的表达图画成一排方便纵向比较fig, axes plt.subplots(2, 5, figsize(20, 8)) for i, (name, ad) in enumerate(adata_list.items()): ax axes[i // 5][i % 5] sc.pl.spatial(ad, colorB2M, axax, showFalse, titlename) ax.set_title(f{name} - B2M, fontsize10) plt.tight_layout() plt.show()如果你做的是同一个组织连续切片的比较建议严格统一vmin和vmax否则每张图各按自己的最大最小值上色会被误读成表达了明显差异。5. 预处理的坑我替你们踩过这些十张切片实战中的血泪总结这部分不是教科书内容全是实操撞出来的。按条列清楚遇到类似问题直接对号入座。5.1 线粒体比例阈值设得太严健康区域被误杀有一张切片坏死区域确实多我就把pct_counts_mt阈值设到 10%。结果过滤完之后整张切片只剩 40% 的 spot而且肿瘤富集区域的信号明显变少。后来一查发现这个组织的某些区域透化程度偏高本底就是 15% 左右但那里的基因表达是有意义的生物信号。调整策略阈值放宽到 25%并且优先观察过滤前后聚类是否保持稳定而不是死磕一个通用阈值。经验是线粒体比例阈值的“最优值”永远取决于你的组织和切片质量不是取决于网上的教程。两张切片分开设阈值优于十张切片统一设一个阈值。5.2 内存爆炸十张切片合并之后 AnnData 不响应了十张 Visium 切片大概多少个 spot以我的数据来说每张 3000-5000 个 spot合并后大约 4 万个 spot × 2000 个高变基因。矩阵本身不算夸张但如果你把全部表达矩阵包括非高变基因都留在内存里再跑聚类和可视化就会非常卡。我的处理是高变基因筛选后直接用adata_all adata_all[:, adata_all.var[highly_variable]]把非高变基因丢掉矩阵小了接近 60%。如果还想更快把表达矩阵转为稀疏格式存储。5.3 坐标比例失调为什么 spot 分布图看起来像被拉扁了matplotlib 默认的 axis 比例是等比的但如果你在一张很大的画布上画很多子图aspect设置不对会出现坐标被拉伸的情况组织看起来要么变扁要么被拉长。scanpy.pl.spatial内部对单图处理得没问题但用subplot时坐标轴比例可能漂移。解决方法是添加一行ax.set_aspect(equal)或者直接在sc.pl.spatial里传axax后手动设置。别小看这一行审稿人看图会非常在意组织形态是否变形。5.4 聚类结果和形态学位置对不上问题多半出在归一化而非聚类有一次聚类出来两个 clusters 在空间图上完全交融没有解剖学规律性。第一反应是聚类分辨率太低调高一圈还是不行。回头看才发现切片合并后我用target_sum1e4做了全局归一化但其中一张切片测序深度特别浅只有别的切片的一半相当于这张切片所有 spot 的表达谱被“稀释”了归一化后噪声放大导致这张切片的聚类信号压不过技术噪声。后来按切片分别归一化的数据重新合并聚类结果立刻变得干干净净。这个坑提醒我数据质控时一定要检查切片间的总 UMI 分布不要在差异可控之前就合并分析。6. 写在最后这套流程用顺手之后还能怎么扩展十张切片的预处理流程沉淀下来之后后续的扩展方向大概有三个一是把空间邻域特征比如每个 spot 周围 neighbor 的平均表达加入分析用squidpy的邻域富集分析快速找到特定细胞型在空间上的聚集模式二是做空间可变基因筛选找出“在组织空间上存在梯度表达”的基因这类基因往往是区域功能差异的新候选三是把预处理和可视化流程封装成一个可复用脚本下次来任何新切片改一下路径就能直接出 QC 报告和空间图。我自己的体会是任何一套流程的价值都不是跑通一次而是能让你下次处理类似数据的时候把精力放在生物学问题上而不是困在坐标缩放和阈值调试里。这套基于 Python 的流程目前在我们手里已经处理过不止一类组织稳定性我比较有底气。

十张10x Visium切片的空间转录组Python预处理与可视化实战

相关推荐

GB28181码流收发实战：RTP/PS封装解析与排障技巧

泡泡地图+AI：快速搭建整本书大纲的几何方法

Redis主从复制“三次握手”源码拆解：状态机与同步机制

更多资讯

基座模型重塑自动驾驶：架构拆解与工程落地实践

Hadoop完全分布式集群从零搭建：配置解析与踩坑排查指南

流批一体不只有Flink，实时数据模型才是真正的分水岭

基于Python与CNN卷积神经网络的会飞昆虫识别

开源数据同步中间件实战：CDC增量、断点续传与多源接入指南

跌倒检测AI系统实战：从六类行为数据集到边缘部署全流程

今日更新

电商智能运营体系从0到1：架构设计、落地步骤与避坑指南

压缩空气储能与热网热惯性联合优化调度：热电联产系统建模与求解

若依框架AI代码生成不稳定？用Agent Skill三层结构+五步法搞定

本周精选

金融AI落地指南：从调研报告到技术选型与工程实践

大模型上下文成本与注意力优化技术深度解析——为何O(n²)平方复杂度无法被常规优化消除｜TaoToken 统一 Key 通道实测

基于GMM的语音识别系统原理与MATLAB实现详解

本月热门

我发现了一个新思路：用 Remotion + Claude Code 像写代码一样自动化生成短视频

Windows下 Codex 中 Chrome 和 Computer Use 插件不可用问题排查及解决参考方式：TaoToken 统一 Key 配置与验证

2026 大模型集体涨价：用 Python 做企业 Token 成本测算与选型避坑（附配置）

读完文章，想聊聊您的网站？