如何基于空间转录组分割图像与基因XY坐标构建基因表达矩阵

核心逻辑 #

你通过cv2.connectedComponentsWithStats得到的output矩阵本质是分割掩码：每个像素的值对应所属的分割对象ID，其中0固定为背景，1~nb_components-1为你要的细胞/对象ID。只需要将每个spot的坐标映射到掩码上拿到对应的对象ID，再聚合生成表达矩阵即可。

分步实现代码 #

1. 导入依赖库 #

import cv2
import numpy as np
import pandas as pd
from scipy.sparse import csr_matrix

2. 定义坐标匹配函数 #

def get_cell_id(x, y, seg_mask):
    # 亚像素浮点坐标转整数像素坐标
    x = int(round(x))
    y = int(round(y))
    # 校验坐标是否在分割图范围内
    if 0 <= y < seg_mask.shape[0] and 0 <= x < seg_mask.shape[1]:
        cell_id = seg_mask[y, x]
        # 过滤背景ID 0
        return cell_id if cell_id != 0 else np.nan
    return np.nan

3. 批量匹配所有spot到对应细胞 #

# 读取你的spot表格，替换为你的实际文件路径
df_spots = pd.read_csv("spots_info.csv")
# output是你之前调用cv2.connectedComponentsWithStats返回的分割掩码矩阵
df_spots["cell_id"] = df_spots.apply(
    lambda row: get_cell_id(row["x"], row["y"], output),
    axis=1
)
# 过滤没有匹配到细胞的spot
df_matched = df_spots.dropna(subset=["cell_id"]).copy()
df_matched["cell_id"] = df_matched["cell_id"].astype(int)

4. 生成表达矩阵 #

普通稠密矩阵（pandas DataFrame格式，缺失值为NaN）

expr_matrix = df_matched.pivot_table(
    index="target_name",  # 行是基因名
    columns="cell_id",    # 列是细胞ID
    values="intensity",   # 填充值为spot强度
    aggfunc="sum"         # 同一细胞内同一基因的多个spot默认求和，可按需改为mean、count等
).fillna(np.nan)

稀疏矩阵（适合大规模数据，节省内存）

# 先将NaN替换为0，再生成CSR格式稀疏矩阵
sparse_matrix = csr_matrix(expr_matrix.fillna(0).values)
# 单独保存行名（基因名）和列名（细胞ID）
gene_list = expr_matrix.index.tolist()
cell_list = expr_matrix.columns.tolist()

注意事项 #

• 坐标一致性校验：如果你的spot坐标输出规则是x对应行、y对应列，需要把get_cell_id里的seg_mask[y, x]调整为seg_mask[x, y]，避免匹配错位
• 聚合规则可自定义：如果需求是统计每个细胞内每个基因的spot数量，把values改为Spot，aggfunc改为count即可
• 大表优化：如果spot数量超过10万，用apply遍历效率较低，可以直接用numpy向量化操作提速：

  x_coords = np.round(df_spots["x"].values).astype(int)
  y_coords = np.round(df_spots["y"].values).astype(int)
  valid_mask = (y_coords >=0) & (y_coords < output.shape[0]) & (x_coords >=0) & (x_coords < output.shape[1])
  cell_ids = np.full(len(df_spots), np.nan)
  cell_ids[valid_mask] = output[y_coords[valid_mask], x_coords[valid_mask]]
  cell_ids[cell_ids == 0] = np.nan
  df_spots["cell_id"] = cell_ids
  

内容的提问来源于stack exchange，提问作者minkprime