相关散点图

# ============================================================
# NK1R+ 血管密度 vs 胶原面积 相关性散点图（模拟数据）
# 风格参考：Cell 2022 FOLR2/CD8A 相关图（ggpubr::ggscatter + stat_cor）
# 输出：PDF（矢量）+ SVG
# ============================================================
rm(list = ls())

# 安装守卫 ---------------------------------------------------
if (!requireNamespace("ggplot2", quietly = TRUE))  install.packages("ggplot2")
if (!requireNamespace("ggpubr",  quietly = TRUE))  install.packages("ggpubr")

library(ggplot2)
library(ggpubr)

set.seed(20261018)

# 1. 模拟数据 -------------------------------------------------
# 设计：胶原面积随纤维化程度（Normal/Mild/Moderate/Severe）递增，
#       NK1R+ 血管密度与胶原面积正相关（r ≈ 0.75），
#       模拟"二者密切相关、且与纤维化程度一致"的结果。
n_per_stage <- 30
fibrosis_stage <- rep(c("Normal", "Mild", "Moderate", "Severe"),
                      each = n_per_stage)

# 纤维化程度驱动的潜在梯度（决定两个变量的共同趋势）
stage_effect <- rep(c(-1.5, -0.5, 0.5, 1.5), each = n_per_stage)

# NK1R+ 血管密度（vessels/mm^2）：受潜在梯度驱动 + 个体噪声
vessel_density <- 40 + 12 * stage_effect +
  rnorm(length(stage_effect), mean = 0, sd = 8)
vessel_density <- pmax(vessel_density, 5)  # 下限截断，保持生物学合理

# 胶原面积（%）：与血管密度共享潜在梯度，故二者强相关
collagen_area <- 25 + 10 * stage_effect +
  0.35 * (vessel_density - mean(vessel_density)) +
  rnorm(length(stage_effect), mean = 0, sd = 6)
collagen_area <- pmin(pmax(collagen_area, 2), 85)  # 限制在 [2, 85] %

data <- data.frame(
  NK1R_vessel_density = round(vessel_density, 1),
  collagen_area = round(collagen_area, 1),
  fibrosis_stage = fibrosis_stage
)
data$fibrosis_stage <- factor(data$fibrosis_stage,
                              levels = c("Normal", "Mild", "Moderate", "Severe"))

# 保存数据（CSV，方便复用）
write.csv(data, "14.1_NK1R_vessel_density_vs_collagen_area.csv",
          row.names = FALSE)

# 快速核对相关系数
cat("Pearson r =", round(cor(data$NK1R_vessel_density,
                             data$collagen_area), 3), "\n")
print(table(data$fibrosis_stage))

# 2. 绘图（Cell 同款：ggscatter + stat_cor + geom_smooth） -----
# 配色：按纤维化程度的连续渐变（与"跟纤维化程度密切相关"呼应）
stage_colors <- c("Normal"   = "#f7fcf5",
                  "Mild"     = "#c7e9c0",
                  "Moderate" = "#74c476",
                  "Severe"   = "#005a32")

p <- ggscatter(data,
               x = "NK1R_vessel_density",
               y = "collagen_area",
               color = "fibrosis_stage",
               palette = unname(stage_colors),
               size = 3,
               add = "none",
               conf.int = FALSE) +
  # 相关系数（pearson，数据已近似线性且无极端离群点）
  stat_cor(method = "pearson",
           label.x = 10, label.y = 52,
           size = 7, label.sep = "\n") +
  # 回归线：黑色实线
  geom_smooth(method = "lm", color = "black", linetype = "solid",
              linewidth = 0.8, se = FALSE) +
  labs(x = "NK1R+ vessel density (vessels/mm²)",
       y = "Collagen area (%)",
       color = "Fibrosis stage") +
  theme_bw() +
  theme(
    # 加粗边框
    panel.border = element_rect(color = "black", linewidth = 0.8, fill = NA),
    # 加粗坐标轴文字/标题
    axis.text  = element_text(face = "bold", size = 14),
    axis.title = element_text(face = "bold", size = 16),
    # 图例
    legend.text  = element_text(face = "bold", size = 11),
    legend.title = element_text(face = "bold", size = 12),
    legend.position = "right",
    # 去掉网格线
    panel.grid.major = element_blank(),
    panel.grid.minor = element_blank()
  )

# 3. 输出 PDF + SVG ------------------------------------------
ggsave(filename = "NK1R_vessel_density_vs_collagen_area.pdf",
       plot = p, width = 5.5, height = 5, device = "pdf")
ggsave(filename = "NK1R_vessel_density_vs_collagen_area.svg",
       plot = p, width = 5.5, height = 5, device = grDevices::svg)

cat("已输出: PDF + SVG\n")