📚 课程目录第2章数据准备

数据准备：理解你的Spot计数矩阵和空间坐标

做批次效应消除，第一步不是调包，而是搞清楚你手里到底有什么数据。我见过不少同学拿到10x Visium的数据就急着跑Harmony，结果跑了半天发现坐标压根没对齐——这种坑踩一次就够你受的。

这一章我们就来手把手搞定数据准备。说白了，就是搞清楚两个东西：计数矩阵和空间坐标。你想想看，没有坐标你连Spot在哪儿都不知道，还谈什么消除批次效应？

两个主流平台：10x Visium vs Slide-seq

目前市面上最常用的两个平台，一个是10x Visium，另一个是Slide-seq。它们的数据结构差异还挺大的，我刚开始切换平台时也翻过车。

特性	10x Visium	Slide-seq
Spot直径	55 μm	10 μm
Spot间距	100 μm	～15 μm
典型Spot数	～5,000	～50,000+
输出格式	H5 + 图像 + JSON	Bead位置CSV + 计数矩阵
坐标系统	组织切片图像坐标(像素)	物理坐标(x, y in μm)

我个人习惯用10x Visium做初筛，因为它分辨率低但组织形态完整；Slide-seq则适合做精细定位，不过数据量大得多——处理50万个Spot的矩阵，内存不够是真的会崩的。

计数矩阵里到底装了啥？

计数矩阵，行是基因，列是Spot。就这么简单。但里面的数值是UMI计数，不是RPKM或者TPM。

  关键点：UMI计数是整数，代表每个Spot里检测到某个基因的独特转录本数。别去归一化，至少在批次效应消除之前别动它。

举个例子，一个10x Visium的标准输出长这样：

         Spot_1 Spot_2 Spot_3
Gene_A     12      0      5
Gene_B      3      8      0
Gene_C      0      1      9

注意看Spot_2的Gene_A是0，这不代表基因没表达——可能是没捕获到。空间转录组的数据稀疏性很高，Spot之间差异巨大，这是常态。

实战经验：我在处理一个肝脏样本时，发现边缘Spot的UMI总数只有中心区域的1/10。这不是批次效应，是组织边缘的捕获效率天然低。后面做批次校正时得小心别把它当批次抹掉了。

空间坐标文件怎么读

10x Visium的坐标存储在tissue_positions_list.csv里。这个文件有6列：

barcode,in_tissue,array_row,array_col,pxl_col_in_fullres,pxl_row_in_fullres
AAACCCAAGAAACAAG-1,1,0,0,543,1212
AAACCCAAGAAACAAT-1,1,1,0,645,1210
AAACCCAAGAAACAAC-1,0,2,0,747,1208

其中in_tissue=1表示这个Spot在组织切片上，in_tissue=0的Spot是背景区域——我在第一次分析时就忘了过滤掉这些背景Spot，导致结果一团糟。

Slide-seq的坐标格式不太一样。它通常给一个beads.csv：

barcode,x,y
AAACCCAAGAAACAAG-1,1200.5,3400.2
AAACCCAAGAAACAAT-1,1215.3,3398.7

这里的x和y是物理坐标，单位通常是μm。Slide-seq没有图像参考，所以坐标就是纯粹的数值位置。

数据读取实操

我习惯用Seurat或Scanpy来读数据。这里给个Python版本的例子，毕竟现在单细胞生态越来越Python化。

读取10x Visium的H5文件：

import scanpy as sc

adata = sc.read_10x_h5("filtered_feature_bc_matrix.h5")
adata.var_names_make_unique()

# 读取空间坐标
import pandas as pd
positions = pd.read_csv("tissue_positions_list.csv",
                         header=None,
                         names=["barcode","in_tissue","array_row",
                                "array_col","pxl_col","pxl_row"])

# 只保留组织上的Spot
positions = positions[positions["in_tissue"] == 1]

# 对齐barcode顺序
positions.index = positions["barcode"]
positions = positions.loc[adata.obs_names]

adata.obsm["spatial"] = positions[["pxl_col","pxl_row"]].values

你会看到最后一步是把坐标放进了obsm里——这是Scanpy的惯例。我建议统一命名为"spatial"，后面很多工具都认这个key。

坑我踩过：千万别用raw_feature_bc_matrix.h5！那个是没过滤的，包含空液滴和背景。一定要用filtered_feature_bc_matrix.h5。我第一回做分析时没注意，多花了两天debug。

理解坐标系的差异

10x Visium的坐标是像素坐标，图片分辨率不同会导致数值范围不一样。Slide-seq给的是物理坐标，跨样本时一般不需要缩放。

这里有个细节你可能会忽视：坐标原点。10x Visium的原点在图像左上角，Slide-seq的原点一般在左下角或者中心。

为什么要关心这个？因为跨批次合并时，坐标必须对齐到同一个参考系。我自己的做法是：

将所有样本的坐标缩放到0-1之间
保留原始坐标用于可视化，批次校正时用归一化坐标
注意检查坐标的宽高比，有些图像是畸形的

样本合并前的检查清单

在合并多个样本之前，我会过一遍这个清单：

所有样本的基因名格式统一（比如都大写）
Spot的barcode格式一致（有些版本带-1有些不带）
坐标列名统一（pxl_col vs x）
没有重复的barcode（不同样本可以有相同barcode，但要加样本前缀）
UMI计数是整数，没有提前归一化

关于barcode重复的问题，我是这样处理的：每个样本加载后立即在obs里加一列sample_id，然后把barcode改成"sample_id-barcode"的形式。这样合并后保证唯一。

结构总览：数据准备的核心步骤

下面这张图总结了从原始数据到可分析的AnnData对象的完整流程。每个步骤之间都有依赖关系，我建议你把它存下来当参考。

关于稀疏矩阵的存储

空间转录组的计数矩阵有多稀疏？我见过一个10x Visium样本，5,000个Spot × 20,000个基因，非零元素不到5%。这种矩阵用常规的密集矩阵存，内存直接爆炸。

解决办法是用稀疏矩阵格式。Scanpy默认就用scipy.sparse.csr_matrix，你不用额外操心。但如果你自己写数据加载，记得指定dtype=np.float32，能省一半内存。

嗯，这里要注意：有些工具在处理稀疏矩阵时会偷偷转成密集矩阵，比如某些聚类算法。我遇到过一次，128GB的机器直接OOM，排查了半天才发现是这个问题。

做一批样本的经验数据

我最近处理了6个10x Visium的免疫组库样本，整理了一些实际数字供你参考：

指标	典型值	备注
每个样本Spot数	～4,500	过滤后
每个Spot UMI数	1,000–10,000	中位数约3,500
检测到的基因数	～18,000	Ensembl版本影响
合并后的AnnData	～400 MB	稀疏矩阵格式
数据加载时间	30秒–2分钟	取决于磁盘速度

注意：如果你发现个别样本的Spot UMI数整体偏低（比如中位数只有500），别急着做批次校正。先检查切片质量——可能是组织破裂或者捕获芯片出了问题。我有个同事花了两周调参数，最后发现是那个样本的切片切碎了。

数据准备到这一步基本就到位了。你拿到了正确的计数矩阵，坐标也对齐了，样本也合并好了。下一件事就是看看这些数据到底长什么样——可视化检查是必不可少的。