写这篇笔记之前先说清楚这是本系列第三篇前两篇我分别讲了单细胞转录组的基础概念、Python生态里常用的数据结构以及分析流程的整体鸟瞰。这一篇不再重复基础而是沿着一条真正能跑的代码路径把空间转录组spatial transcriptomics分析流程从头到尾拆一遍。我默认你手里已经有一份10x Visium的空间转录组数据无论是公开数据集还是自己的测序结果。如果你还没接触过单细胞建议先回去补前两篇否则后面这些代码和概念叠在一起容易晕。本篇的重点是环境怎么搭、数据怎么读、质控怎么做、标准化聚类差异分析怎么串起来、最后怎么在空间维度上把结果“画”出来。每一步我都会贴实际代码和参数解释顺带分享我踩过的坑。1. 为什么第三篇才真正“跑起来”从单细胞到空间转录组的流程骨架1.1 先分清单细胞与空间转录组的分析差异单细胞转录组和空间转录组在分析思路上有一个根本区别单细胞数据通常只有“细胞身份”你拿到的是一个细胞一个表达谱空间转录组则多了一个“坐标”每个spot测序位点不仅有自己的表达谱还有在组织切片上的XY坐标。这个差异直接决定了分析流程的走向。在单细胞流程里我们做完质控、标准化、聚类、注释后任务基本就结束了最多做个细胞通讯、拟时序。但在空间转录组里聚类结果出来后还要回到坐标上去看某个细胞类型在空间上是不是扎堆某两个细胞类型的邻域关系是不是显著富集某些基因的表达是不是沿空间梯度变化这就需要在传统单细胞流程之上额外叠加一套空间分析模块。另外数据分辨率的差异也要提前有心理准备。10x Visium每个spot捕获的是周围几个细胞的混合转录本不是单细胞分辨率。这意味着表达矩阵里每个样本是spot而不是细胞spot数量通常只有数千到几万个而单细胞动辄十万级。后续做QC阈值、聚类参数时都要围绕spot的特性去调整。1.2 全流程模块梳理从原始数据到结果输出把整个流程拆成模块来看跑起来之后其实是有清晰依赖关系的数据读取把10x Visium的原始输出h5或mtx矩阵 空间坐标读成AnnData对象。质控过滤计算每个spot的基因数、UMI数、线粒体基因占比过滤掉质量差的spot和低表达的基因。标准化与对数化消除测序深度差异让spot之间可比。特征选择挑出高变基因降低下游计算量。降维PCA将高维表达压缩到几十个主成分。邻域图构建与聚类基于PCA空间构建k近邻图用leiden/louvain算法划分spot群。差异分析比较不同spot群之间的表达差异提取marker基因。空间分析将聚类结果、基因表达映射到组织坐标上用Squidpy这类工具做空间邻域富集、空间表达模式分析。这个流程和你熟悉的GEO数据挖掘的“下载、处理、质控、差异分析”框架是相通的核心差别在于单细胞/空间平台有自己的数据格式和QC指标后续聚类、注释环节的权重更高。我前两篇一直在强调“流程感”这一篇就是把这个感觉落实到代码里。2. Python环境搭建与库选型一次装对少走三天弯路2.1 用conda建独立环境避免污染空间转录组分析涉及的Python库相当多而且很多库对版本极度敏感比如anndata、scanpy、numpy之间如果版本不对读数据时会直接报一堆莫名其妙的错。因此我强烈建议不要用base环境而是用conda单独建一个环境。conda create -n spatrans python3.9 -y conda activate spatransPython版本我推荐3.9。3.8偏旧部分新库已经不再支持3.10以上某些库尤其是老版本的leidenalg容易出现编译兼容问题。3.9是当前单细胞生态兼容性最好的中间地带我实测下来最稳。建好环境后再装库安装顺序也很重要先装scanpy主轴再装辅助库避免pip自动解析依赖时发生冲突。2.2 必须装齐的核心库与版本不同数据规模、不同分析目标需要的库略有差异但基本盘是这七个scanpy分析流程主框架读写、QC、标准化、聚类、差异分析全靠它。anndata数据容器库scanpy底层依赖存储表达矩阵和元数据。squidpy空间分析专用库邻域富集、空间统计、形态学分析都用它。leidenalg / python-igraphleiden聚类算法的后端scanpy聚类时要用。matplotlib / seaborn绘图基础库。pandas / numpy / scipy数据处理基础设施必装。scikit-learn偶尔用它的层次聚类、PCA等接口作补充。安装命令我一般这样写pip install scanpy anndata squidpy1.4.1 leidenalg python-igraph matplotlib seaborn numpy pandas scipy scikit-learn这里建议把squidpy单独装因为它和scanpy的版本匹配比较敏感。如果直接pip install squidpy它会把scanpy升级到最新版而最新版scanpy有时会引入不兼容改动。稳妥做法是先装scanpy再装squidpy装完用两条命令验证一下版本import scanpy as sc import squidpy as sq print(sc.__version__, sq.__version__)2.3 安装时的常见坑源替换与依赖冲突国内用户最常遇到的是下载慢、超时。我的做法是配置pip国内镜像源实测速度能差出几十倍pip config set global.index-url https://pypi.tuna.tsinghua.edu.cn/simple如果你已经装了一半出问题比如提示“Cannot uninstall ‘numpy’”或者某个依赖出现冲突不要硬着头皮继续pip install建议直接把环境删掉重建省时省力conda deactivate conda remove -n spatrans --all -y conda create -n spatrans python3.9 -y另外提醒一句不要在root用户下直接改全局环境也不要用sudo pip很容易把系统级的Python搞坏。保持conda环境隔离后患少很多。3. 数据读取与质控实操喂给流程的第一口原料3.1 10x Visium数据长什么样10x Visium是目前空间转录组最主流的数据平台它的官方Cell Ranger输出目录一般长这样outs/ ├── filtered_feature_bc_matrix.h5 ├── filtered_feature_bc_matrix/ │ ├── barcodes.tsv.gz │ ├── features.tsv.gz │ └── matrix.mtx.gz ├── spatial/ │ ├── tissue_positions_list.csv │ ├── scalefactors_json.json │ └── tissue_hires_image.pngfiltered_feature_bc_matrix.h5是过滤后的表达矩阵包含了基因名、spot条形码和计数spatial文件夹里tissue_positions_list.csv是每个spot在组织切片上的坐标和是否在组织范围内的标记scalefactors_json.json是低分辨率/高分辨率图像与全分辨率坐标之间的换算比例。整个分析中表达矩阵和空间坐标必须一一对应如果你的数据是自己处理的务必保证barcode名称完全一致这是后续所有空间分析的地基。3.2 read_visium读取代码与参数说明scanpy封装好了Visium数据的读取函数一行代码就能把表达矩阵和空间坐标读进来import scanpy as sc adata sc.read_visium( pathouts/, count_filefiltered_feature_bc_matrix.h5, library_idCytAssist_FFPE_Human_Lung, load_imagesTrue, )这里library_id是什么它对应spatial目录下tissue_positions_list.csv所属的样本标识。如果只有一个样本scanpy会自动识别但你最好显式指定否则后续做多样本合并时容易搞混。load_imagesTrue会读取组织图像后面space图要把spot画在图像上这一步必须开启。读取后顺手检查一下数据结构print(adata.shape) print(adata.obs.head()) print(adata.uns[spatial])adata.shape显示的是(spot数, 基因数)正常情况下spot数应该在几千到数万之间基因数在两万上下。obs里应该有array_row、array_col芯片上的行列号、pxl_row_in_fullres、pxl_col_in_fullres全分辨率像素坐标这些字段在后面画空间图时会被scanpy自动调用。3.3 QC指标计算与阈值选择拿到原始矩阵后第一步永远是质控。这里和单细胞类似但阈值要按spot的特性来调整。我先说两个通用指标的计算adata.var[mt] adata.var_names.str.startswith(MT-) sc.pp.calculate_qc_metrics( adata, qc_vars[mt], percent_topNone, log1pFalse, inplaceTrue, )calculate_qc_metrics会在obs里生成n_genes_by_counts每个spot检测到的基因数、total_counts每个spot的总UMI数、pct_counts_mt线粒体基因UMI占比等列。为什么要算线粒体基因占比这是一个很经典的质量信号如果一个spot的细胞破裂细胞质里的mRNA大量流失线粒体基因相对占比会异常升高。在单细胞中常用的阈值是pct_counts_mt20%但空间转录组spot捕获的是多个细胞的混合线粒体基因占比会相对偏低一般建议先看分布再定经验值可以放在5%-15%之间。过滤操作建议做两轮先粗筛spot再筛基因。sc.pp.filter_cells(adata, min_genes200) sc.pp.filter_genes(adata, min_cells3)min_genes200的意思是检测到基因数少于200的spot直接去掉。与单细胞不同每个spot的建库效率通常低于单细胞基因数少于200的spot大概率是芯片边缘或者组织碎裂区域。min_cells3表示只在少于3个spot里表达的基因丢掉这类基因要么是背景噪音要么是极低表达留着只会增加后面多重检验校正的负担。阈值怎么科学地定我习惯先画分布图再看数据形态。比如把n_genes_by_counts画成小提琴图、把pct_counts_mt画成散点图根据分布用中位数加减3倍MAD绝对中位差来定阈值。快速看一眼import matplotlib.pyplot as plt sc.pl.violin(adata, keys[n_genes_by_counts, total_counts, pct_counts_mt], multi_panelTrue) plt.show()如果小提琴图显示大部分spot的基因数集中在150-200之间那么min_genes200可能过滤太多需要下调到100如果集中在300以上200就是合理的。3.4 过滤与可视化检查过滤一定要链式写避免中间变量覆盖adata.obs[n_genes_by_counts] adata.obs[n_genes_by_counts].astype(int) adata.obs[total_counts] adata.obs[total_counts].astype(int) adata adata[adata.obs[n_genes_by_counts] 200, :] adata adata[adata.obs[pct_counts_mt] 15, :]过滤完再看一眼形状心里要有数print(f过滤后剩余spots{adata.n_obs}) print(f过滤后剩余genes{adata.n_vars})过滤本身并不难难的是确认过滤没有把不该丢的数据丢掉。我习惯过滤前用sc.pl.spatial画一次总UMI数的空间分布看看低质量的spot是不是集中在组织边缘。如果组织内部也有一片低质量spot那可能是切片质量或芯片问题而不是正常的背景区域这时候盲目过滤会丢失真正的组织信号。质控不是机械地套阈值而是要看空间分布。4. 标准化、高变基因、降维、聚类与差异分析流程中最核心的五步4.1 标准化为什么默认 normalize_total log1p数据过滤之后叠加上TCGA/GEO数据挖掘的做法会更直观芯片数据要做背景校正、归一化单细胞/空间转录组也类似但要处理的批次效应来自测序深度差异。不同spot捕获到的RNA总量不同总UMI数可能差出好几倍如果不做标准化聚类结果会被测序深度主导真实组织结构容易被淹没。scanpy两步走是教科书级别的标配sc.pp.normalize_total(adata, target_sum1e4) sc.pp.log1p(adata)normalize_total做的事情是把每个spot的counts按比例缩放让每个spot的总UMI数都等于target_sum通常取1e4。为什么要选1e4这是一个经验约定相当于把每个spot的测序深度拉到同一个量级后续表达值近似于“每万个UMI中该基因有几个计数”。log1p则是取log(x1)一方面压缩动态范围把表达量从0-几万压缩到0-10左右的区间另一方面让低表达基因的差异不再被高表达基因压制。这里必须强调一个常见错误你的adata.X必须是原始整数计数才能做normalize_total绝不能对已经log过的矩阵再做一次normalize和log1p。如果数据是别人处理好的log1p矩阵你拿到手后应该直接从PCA开始而不是重复标准化。判断方法很简单如果adata.X里的最大值超过20甚至上千大概率是原始计数如果最大值在10左右且有0-1之间的浮动大概率已经log过。拿不准时可以看raw部分如果在读入数据时设了adata.raw adata.copy()那就用raw里的原始X来做标准化。4.2 高变基因与PCA把“噪音”降下来标准化后不要直接聚类矩阵里有两万多个基因其中大部分在spot之间波动不大聚类时反而是噪音。所以要挑高变基因。scanpy提供了现成接口sc.pp.highly_variable_genes(adata, n_top_genes2000, flavorseurat_v3, layercounts)flavorseurat_v3需要把原始counts存在layer里它用的是基于方差均值比的方法。如果你更熟悉scanpy默认的flavorseurat也可以直接用效果差异不算太大。n_top_genes通常取2000到3000取太少会漏掉一些弱信号取太多又引入噪音2000是大多数分析里比较好的默认值。把高变基因筛选出来之后adata adata[:, adata.var[highly_variable]].copy() sc.tl.pca(adata, svd_solverarpack, n_comps50)PCA是一种线性降维把几千个高变基因的表达模式压缩成50个左右的主成分每个主成分都是原基因的线性组合。这个步骤的意义是后续计算spot之间的相似度时不再用几千维的表达向量而是用几十维的主成分得分计算量大幅下降同时还能去噪。svd_solverarpack适合矩阵较大但不需要计算全部主成分的场景n_comps50也可以用默认值scanpy会根据数据规模自动调整。降维后画一张碎石图看前50个主成分里还保留了多大比例的方差。如果前10个主成分已经解释了70%以上的方差后续neighbors里n_pcs取20就够了如果解释率偏低则可以考虑保留30-40个主成分。4.3 邻域图与leiden聚类resolution怎么定PCA之后是构建邻域图和聚类这是把spot分群的关键环节sc.pp.neighbors(adata, n_neighbors15, n_pcs30) sc.tl.leiden(adata, resolution0.8, flavorigraph, n_iterations2, directedFalse)neighbors做的事是在PCA降维后的空间里为每个spot找到15个相似度最高的邻居构建一个加权k近邻图。n_neighbors越小图的连通性越弱聚类得到的群越碎越大则越容易合并成大群。n_pcs30表示用前30个主成分来计算距离如果前面PCA的方差解释率不够可以适当调高到40。leiden聚类是在邻域图上做社区发现。leiden算法是对louvain的改进能保证输出分区在连接性上有数学保证。它的核心参数是resolution值越大社区划分越细碎群数越多。实际操作经验0.5-1.0之间是探索性分析常用的区间。我一般会跑一组resolution的网格观察聚类群数从5到20的变化然后根据是否匹配组织学结构来判断。for res in [0.2, 0.5, 0.8, 1.2]: sc.tl.leiden(adata, resolutionres, key_addedfleiden_{res})不要把希望全部寄托在默认参数上。科研项目里聚类数不是算法“算出来”的而是你要根据marker基因和组织学先验知识去选择一个生物学上合理的尺度。resolution0.8跑出一个leiden列不代表这是唯一答案。4.4 层次聚类辅助判断聚类数的一个补充手段层次聚类在单细胞/空间转录组里经常被忽视但我认为它至少有两个价值一是可以在选择resolution时快速评估数据自然的类簇数量二是当你对某几个spot群的关系存疑时可以画一张树状图看它们的层次关系。scanpy本身没有直接封装层次聚类到AnnData的接口但借助scipy很容易实现。把PCA降维后的表达矩阵拿出来用层次聚类算法做聚合再用树状图辅助判断import numpy as np from scipy.cluster.hierarchy import linkage, dendrogram, fcluster from scipy.spatial.distance import pdist X_pca adata.obsm[X_pca][:, :30] dist_mat pdist(X_pca, metriceuclidean) link_mat linkage(dist_mat, methodward, optimal_orderingTrue) dendrogram(link_mat, no_labelsTrue, color_threshold0.7 * max(link_mat[:, 2])) plt.show()ward法是一种通过最小化合并后类内方差的方式来进行层次聚类的方法它在单细胞分析中通常表现较好因为对噪声相对鲁棒。看到树状图后你可以根据最长的分支数目大致判断数据中明显的类簇个数。如果想直接得到离散标签可以用fclusterspatial_labels fcluster(link_mat, t12, criterionmaxclust)这个方法不需要先算neighbors而是直接在PCA距离上用凝聚聚类来做两者结果不会完全一致但如果在某个分辨率下leiden结果和层次聚类的类簇数差异非常大你就要警惕是不是PCA的主成分选择或者neighbors参数设得不合适。4.5 差异分析从聚类结果中提取marker基因聚类完成后下一步就是差异分析这一步和GEO数据挖掘的差异分析逻辑是一样的只是分组变量从“疾病/对照”变成了“leiden聚类群”。scanpy提供了一套很方便的管道sc.tl.rank_genes_groups(adata, groupbyleiden, methodwilcoxon, use_rawTrue)method有几种选择wilcoxonWilcoxon秩和检验是单细胞/空间转录组中最常用的因为它不要求数据符合正态分布对离群值不敏感t-test假设正态分布检验功效更高但更容易受极端值干扰logreg是逻辑回归方法速度慢但结果更稳定。我通常首选wilcoxon。use_rawTrue的意思是在原始的、未做标准化之前的counts矩阵上做检验。为什么因为差异检验应该基于原始计数用log1p后的矩阵做检验虽然也能跑但会引入标准化带来的伪差异。确保读取数据后保存过adata.rawadata.raw adata.copy()执行完rank_genes_groups后结果是乱的要提成DataFramedf_de sc.get.rank_genes_groups_df(adata, group0).sort_values(pvals_adj) print(df_de.head())只看pvals_adj小于0.05的基因并不够实际分析中还要看logfoldchanges的绝对值。单细胞和空间转录组数据量庞大检验功效很高很多基因的差异p值非常显著但实际变化倍数只有0.1这种基因生物学意义不大。我的习惯是先用logfoldchanges1过滤pvals_adj0.05只是第一道门槛markers df_de[(df_de[pvals_adj] 0.05) (df_de[logfoldchanges].abs() 1)]拿到marker基因之后可以用点图或热图快速验证聚类群的生物学身份marker_genes [EPCAM, PTPRC, VWF, COL1A1, CD3D, CD79A, LYZ] sc.pl.dotplot(adata, var_namesmarker_genes, groupbyleiden, use_rawTrue)这一步是连接“数学聚类”和“生物学注释”的桥梁。如果某个群高表达EPCAM基本可以判定是上皮细胞群高表达CD3D/CD79A则是T/B淋巴细胞VWF是内皮LYZ是髓系。注释的结果直接决定了后续所有空间分析的故事走向宁可多花时间读文献验证marker也不要凭一个基因就拍板。5. 空间维度的进阶操作Squidpy与空间可视化5.1 空间邻域构建与富集分析传统单细胞流程做到marker注释基本就算收工了空间转录组还多一道工序把分析结果放回空间坐标挖掘组织微环境信息。Squidpy是这个环节最常用的库它的核心函数可以用两行代码完成空间邻域的构建import squidpy as sq sq.gr.spatial_neighbors(adata, coord_typegeneric, n_neighs6)spatial_neighbors会基于每个spot的物理坐标构建空间图n_neighs6表示只连接欧氏距离最近的6个邻居。Visium芯片上6个邻居刚好对应一个六边形周围的6个spot这个参数在网格型空间数据里非常经典。有了空间邻域图就能做邻域富集分析neighborhood enrichment。它背后的数学逻辑是随机打乱spot群标签计算每种细胞类型对的期望共现频率再和实际共现频率比较得到富集分数sq.gr.nhood_enrichment(adata, cluster_keyleiden) sq.pl.nhood_enrichment(adata, cluster_keyleiden)输出是一个矩阵行和列都是leiden群数值大于0代表两个群倾向于彼此邻近小于0代表空间分离。这张矩阵能告诉你很多生物学信息比如肿瘤区域的上皮细胞和成纤维细胞是否共定址免疫细胞是否集中在某个特定区域。另一个常用功能是共现分析co-occurrence计算某种细胞类型在给定空间尺度下的共现趋势sq.gr.co_occurrence(adata, cluster_keyleiden) sq.pl.co_occurrence(adata, cluster_keyleiden, clusters3)共现分析可以帮助你定量描述“某种细胞随着与另一类细胞的物理距离增加共现概率如何衰减”。这个指标在肿瘤微环境研究里特别有价值。5.2 画空间图的关键参数与踩坑空间转录组分析的颜值担当就是spatial图但scanpy默认参数在实际使用中有很多可以优化的地方。先看最简单的画法sc.pl.spatial( adata, colorleiden, spot_size55, palettetab20, frameonFalse, saveleiden_spatial.pdf, )spot_size控制每个spot在图像上的物理大小数值太大会盖住组织形态细节太小则看不清分布规律。我通常先用55-70试一张再根据图像反馈微调。palette建议指定scanpy默认的色板在群数超过12个时颜色区分度会明显下降tab20是20类颜色的好选择如果群数超过20就考虑用matplotlib的hsv循环色图搭配自定义调色板。frameonFalse可以把图像周围的边框去掉导出后方便直接放在论文/汇报ppt里。save如果不指定则直接在交互窗口弹出指定后自动保存到当前目录的figures文件夹下。还有一个高频坑画图时横坐标太密集。热图或小提琴图的x轴如果样本/变量太多标签会挤成一团看起来像是黑色的一坨这是网友常问的“python画图横坐标太密集”问题。解法有几个import matplotlib.ticker as ticker fig, ax plt.subplots(figsize(8, 4)) sc.pl.dotplot(adata, var_namesmarker_genes, groupbyleiden, axax) ax.xaxis.set_major_locator(ticker.MaxNLocator(nbins6)) plt.xticks(rotation45, haright)本质上就是三件事把画布拉宽、限制显示刻度数量、旋转标签。5.3 单细胞与空间数据联合分析的常见思路手头同时有单细胞数据和空间转录组数据时常见的做法是先用单细胞数据做精细注释因为单细胞分辨率高能识别出稀有细胞类型再把注释结果映射到空间数据上。具体操作是把两份数据的公共基因提取出来用scVI、CellTypist这类工具做标签迁移或用Seurat的label transfer思路在Python里复现先建参考数据的PCA模型再用空间表达矩阵去预测每个spot的细胞类型比例。这本质上是一个监督学习问题。Scanpy本身不直接提供参考映射接口我一般会用scvi-tools来做import scvi scvi.model.SCVI.setup_anndata(adata_sc, batch_keysample) vae scvi.model.SCVI(adata_sc) vae.train()跑完后可以用scvi的“deconvolution”模型估计每个spot的细胞类型混合比例。对于Visium这类低分辨率平台spot里混着多种细胞是常态不做反卷积、直接把整群标记成一个细胞类型容易出大偏差。但这个方法门槛稍高本篇先点到为止下一篇我会专门展开写。6. 常见问题与排查技巧实录实战速查表6.1 六类高频报错与解决方案我在跑流程的半年时间里碰到的问题有大半集中在下面这几个地方做成表格方便你对症下药报错/现象常见原因解决办法ValueError: X has different shapeanndata版本升级后shape不匹配或old expression矩阵和新obs/var对不上检查anndata版本确认所有过滤操作都用adata adata[...]重新赋值避免在子视图上继续改MemoryError内存不足表达矩阵以float64全量载入或spot数太多读入数据后立刻转float32用sc.pp.filter_genes筛掉低表达基因或在读取阶段设置dtypenp.float32leiden聚类结果每次跑都不一样leiden算法有随机性默认没有固定随机种子在sc.tl.leiden中加random_state0并用flavorigraph提升稳定性ImportError: libgtk-3.so.0Linux无图形环境下matplotlib的GTK后端缺失代码开头设置import matplotlib; matplotlib.use(Agg)避免调用GUI后端KeyError: leiden聚类步骤还没跑就尝试画图或key名称拼写不一致检查adata.obs.columns确认聚类结果列名用key_added参数统一取名FileNotFoundError读取Visium数据失败spatial文件夹路径不对或library_id与数据不匹配确认数据目录包含spatial/tissue_positions_list.csv显式指定library_id不要依赖自动识别这里特别提醒第三行关于随机性不只是leiden任何聚类算法在多次运行之间的细微差异都是正常的。你在复现别人代码时如果发现聚类数相同但某些spot的群标签对不上先检查随机种子如果还不行再检查输入矩阵是否被修改过。6.2 可视化与内存问题的独家小技巧针对网友问烂的“python画图横坐标太密集”问题我补一个更精简的通用方案。如果用的是matplotlib最简单的方式是plt.xticks(rotation45, haright) plt.tight_layout()haright是水平对齐方式能让标签旋转后仍然对齐刻度位置tight_layout()自动调整子图边距避免标签被截断。如果你在Jupyter里画完发现横坐标还是挤在一起优先把figsize加大而不是修改字体大小改字体大小往往治标不治本。内存优化方面有一个我实测很有效的小技巧在读取Visium数据之后如果把X转换成了float32再配合scipy的稀疏矩阵存储内存占用量大约能降一半多。Scanpy在默认情况下对于计数矩阵使用CSR稀疏格式比较省内存但PCA之后生成的obsm通常是稠密矩阵spot数乘以主成分数在spot上万之后这个矩阵也会开始吃内存。如果内存告急可以适当降低n_comps比如从50降到30。另外如果你要反复调试参数建议在每次迭代前把过滤后、标准化前的AnnData用zarr格式保存一份之后重新跑流程时直接从这一步加载能省下大量重复读取和质控的时间。代码很简单adata.write_zarr(adata_multi_sample_filtered.zarr)然后下次直接adata sc.read_zarr(adata_multi_sample_filtered.zarr)空间转录组数据动辄几个G多存几个中间版本不算浪费有时候能救回一整天的工作量。结尾写到这儿一条完整的“scanpy squidpy”空间转录组分析主链路已经走通了从数据读取、QC、标准化、PCA、聚类到差异分析最后再叠加空间邻域富集和可视化。我个人最大的体会是这套流程真正难的不是某个函数怎么调而是每一步的“判断”怎么下阈值怎么定、分辨率怎么选、marker怎么验证。代码是死的数据是活的。如果你正在跑自己的数据卡住的时候别硬刚先把四张图同时打出来——小提琴图看质控、UMAP看聚类、spatial图看空间结构、marker点图看注释——四张图放在一起基本能定位80%的问题。最后再分享一个我一直保留的习惯每次修改参数都顺手把关键设置记在代码顶部的注释里包括脚本运行日期、输入数据版本、Python库版本。数据一多版本一乱这些记录就是救命稻草。下一篇我准备专门写多样本整合和批次效应校正那是空间转录组从“能跑”到“跑好”的一道分水岭我们到时候见。