383_单细胞空间重建方法¶

library(spacexr)     # 包含RCTD工具（原名 RCTD 包）
library(Seurat)      # 空间数据处理
library(ggplot2)

# ============================================================
# 步骤1：加载参考 scRNA-seq 数据（带细胞类型注释）
# ============================================================
# 参考数据：已注释好细胞类型的scRNA-seq
ref_seurat <- readRDS("reference_scRNA.rds")

# 提取参考所需信息
ref_counts <- GetAssayData(ref_seurat, layer = "counts")   # 原始count矩阵
ref_celltypes <- factor(ref_seurat$cell_type)               # 细胞类型标签
names(ref_celltypes) <- colnames(ref_seurat)                # 确保名称匹配

# 创建 RCTD 参考对象
ref_object <- Reference(
  counts = ref_counts,       # 基因×细胞 count矩阵
  cell_types = ref_celltypes, # 细胞类型向量
  nUMI = colSums(ref_counts)  # 每个细胞的总UMI数
)

# ============================================================
# 步骤2：加载空间转录组数据（Visium）
# ============================================================
visium_obj <- Load10X_Spatial(
  data.dir = "visium_output/",   # CellRanger Space输出目录
  filename = "filtered_feature_bc_matrix.h5"
)

# 提取空间数据
spatial_counts <- GetAssayData(visium_obj, layer = "counts")
spatial_coords <- GetTissueCoordinates(visium_obj)   # 每个spot的坐标(x,y)
spatial_nUMI   <- colSums(spatial_counts)             # 每个spot的总reads

# 创建 RCTD 空间对象
query_object <- SpatialRNA(
  coords = spatial_coords[, c("imagerow", "imagecol")],  # 空间坐标
  counts = spatial_counts,                                # 表达矩阵
  nUMI   = spatial_nUMI                                  # 总reads
)

# ============================================================
# 步骤3：运行 RCTD 反卷积
# ============================================================
rctd_obj <- create.RCTD(
  spatialRNA = query_object,    # 空间数据
  reference  = ref_object,      # 参考scRNA-seq
  max_cores  = 8                # 8线程并行
)

# 运行反卷积（mode: "multi"=允许多种细胞类型共存，"doublet"=假设每spot最多2种）
rctd_obj <- run.RCTD(rctd_obj, doublet_mode = "multi")

# ============================================================
# 步骤4：提取和可视化结果
# ============================================================
# 提取细胞类型比例矩阵（spot × 细胞类型）
weights_matrix <- as.data.frame(rctd_obj@results$weights)
weights_norm   <- sweep(weights_matrix, 1, rowSums(weights_matrix), "/")  # 归一化

# 将结果整合回 Seurat 对象
for (ct in colnames(weights_norm)) {
  visium_obj@meta.data[[paste0("RCTD_", ct)]] <- weights_norm[, ct]
}

# 空间可视化某种细胞类型的分布
SpatialFeaturePlot(
  visium_obj,
  features = "RCTD_T_cell",     # T细胞的空间分布
  pt.size.factor = 1.5,
  alpha = c(0.1, 1)             # 透明度随值变化
) + scale_fill_viridis_c(name = "T细胞比例")

import tangram as tg    # pip install tangram-sc
import scanpy as sc
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

# ============================================================
# 步骤1：加载数据
# ============================================================
# 参考 scRNA-seq（带细胞类型注释）
sc_adata = sc.read_h5ad("reference_scRNA.h5ad")

# 空间转录组（Visium，10x格式）
sp_adata = sc.read_visium("visium_output/")
sp_adata.var_names_make_unique()
# 标准化空间数据
sc.pp.normalize_total(sp_adata)
sc.pp.log1p(sp_adata)

# ============================================================
# 步骤2：找共同高变基因（减少计算量）
# ============================================================
# 只用两个数据集共有的高变基因做映射
sc.pp.highly_variable_genes(sc_adata, flavor="seurat_v3", n_top_genes=3000)
hvg = sc_adata.var_names[sc_adata.var["highly_variable"]]
common_genes = hvg.intersection(sp_adata.var_names)

# 过滤到公共基因
ad_sc_sub = sc_adata[:, common_genes].copy()   # scRNA-seq子集
ad_sp_sub = sp_adata[:, common_genes].copy()   # 空间子集

print(f"公共高变基因数: {len(common_genes)}")

# ============================================================
# 步骤3：Tangram 映射
# ============================================================
# 预处理（Tangram需要）
tg.pp_adatas(ad_sc_sub, ad_sp_sub, genes=list(common_genes))

# 训练映射模型（优化细胞→spot的映射矩阵）
ad_map = tg.map_cells_to_space(
    adata_sc=ad_sc_sub,      # scRNA-seq参考
    adata_sp=ad_sp_sub,      # 空间数据
    mode="cells",            # 细胞级别映射（也可用"clusters"）
    density_prior="rna_count_based",  # 密度先验（基于spot reads数）
    num_epochs=500,          # 训练迭代次数
    learning_rate=0.1,
    device="cpu"             # 有GPU时用"cuda"更快
)

# ad_map 包含映射矩阵：每个细胞在每个spot的概率

# ============================================================
# 步骤4：投影细胞类型到空间
# ============================================================
# 将 scRNA-seq 的细胞类型信息投影到空间
tg.project_cell_annotations(
    ad_map, ad_sp_sub,
    annotation="cell_type"    # scRNA-seq中的细胞类型列名
)

# 空间可视化（各细胞类型在组织中的分布）
sc.pl.spatial(
    ad_sp_sub,
    color=["T_cell", "B_cell", "Macrophage"],  # 展示多种细胞类型
    size=1.5,
    cmap="Reds",
    save="_tangram_mapping.pdf"
)