# ============================================================
# 04_annotation.R —— 细胞类型注释（对应教程第 09 章）
# 用法：Rscript code/04_annotation.R（或 RStudio 逐段运行）
# 前置：已运行 03_preprocess.R，读取 data/processed/visium_preprocessed.rds
# 说明：本脚本用"marker 基因打分"策略做快速注释；
#       更精细的 cell2location 去卷积见 code/04b_cell2location_squidpy.py（Python）。
# ============================================================

library(Seurat)
library(patchwork)

obj <- readRDS("data/processed/visium_preprocessed.rds")

## 1. 差异表达找每群 marker（辅助判断 cluster 身份）
markers <- FindAllMarkers(obj, only.pos = TRUE, min.pct = 0.25,
                          logfc.threshold = 0.5, assay = "SCT")
top_markers <- markers %>%
  dplyr::group_by(cluster) %>%
  dplyr::slice_max(n = 5, order_by = avg_log2FC)
print(as.data.frame(top_markers))

## 2. 定义细胞类型 marker 集（与教程第 09 章一致）
marker_list <- list(
  "肿瘤细胞"  = c("EPCAM", "KRT8", "KRT19"),
  "CAF"       = c("COL1A1", "COL3A1", "DCN", "PDGFRA"),
  "CD8 T"     = c("CD3D", "CD8A", "GZMB"),
  "巨噬细胞"  = c("CD68", "CD163", "C1QA"),
  "B 细胞"    = c("MS4A1", "CD79A"),
  "内皮细胞"  = c("PECAM1", "VWF")
)

## 3. 模块打分：每个 spot 得到 6 个细胞类型得分
obj <- AddModuleScore(obj, features = marker_list, name = "Score_")
score_cols <- paste0("Score_", seq_along(marker_list))
# 取最高分类型作为该 spot 的"主要细胞类型"
obj$celltype <- colnames(obj[[]])[score_cols][max.col(obj[[score_cols]])]

## 4. 看结果
p1 <- SpatialDimPlot(obj, group.by = "celltype", label = FALSE) +
      ggtitle("主要细胞类型（空间）")
p2 <- DimPlot(obj, group.by = "celltype", label = TRUE) +
      ggtitle("主要细胞类型（UMAP）")
print(p1 | p2)

# 细胞类型构成条形图
comp <- table(obj$celltype)
barplot(comp, las = 2, col = RColorBrewer::brewer.pal(6, "Set2"),
        main = "细胞类型构成", ylab = "spot 数")

## 5. 保存（含 celltype 注释）
saveRDS(obj, "data/processed/visium_annotated.rds")
cat("注释完成，已保存 data/processed/visium_annotated.rds\n")
