# 复现: Nature Communications 2026, doi:10.1038/s41467-026-77288-4, Fig. 4a
# 唾液 sfRNA 抑制 IFN 应答: Huh7 细胞 IFN-β 相对基因表达量(24h/48h, L/M/H 唾液 sfRNA 浓度组)
#
# 数据: 真实数据, 来自论文 Source Data 文件 "Source data for graph.xlsx" 的 'Fig. 4a' sheet
#       (IFN-β relative expression in Huh7 cells infected with saliva; 每组 2 个生物学重复)。
#       图中水平虚线为 mock 感染组值, Source Data 未给出数值, 按原图刻度估计约为 1.0(mock_value=1.0 明确标注为估计值)。
# 显著性字母 a/a/a/ab/b/c 为原文 Fig. 4a 标注的 Fisher LSD 事后检验结果分组(按原图转录, 非本脚本计算)。
# 48h 柱的黑色斜线纹理: ggplot2 无内置纹理填充, 本脚本用 hatch_segments() 手工计算 45° 斜线段绘制。
# 沙箱实测: ~/miniforge3/bin/Rscript reproduce_fig.R -> reproduce_r.png, 无报错。
try(invisible(Sys.setlocale("LC_ALL", "C.utf8")), silent = TRUE)

suppressPackageStartupMessages({
  library(ggplot2)
  library(ragg)
  library(systemfonts)
})

here <- getwd()  # 在 day-2026-10-18-figure 目录下运行本脚本

systemfonts::register_font("noto",
  plain = file.path(here, "NotoSansSC-Regular.ttf"))

df <- read.csv(file.path(here, "data.csv"), stringsAsFactors = FALSE)

COLORS <- c(L = "#EECCF9", M = "#B856D7", H = "#76069A")  # 原图逐像素采样
LETTERS <- c("24h.L" = "a", "24h.M" = "a", "24h.H" = "a",
             "48h.L" = "ab", "48h.M" = "b", "48h.H" = "c")
MOCK_VALUE <- 1.0  # 估计值: 原图虚线位置经刻度校准约 1.0, Source Data 未提供

agg <- aggregate(value ~ timepoint + sfRNA_level, data = df,
                 FUN = function(v) c(mean = mean(v), sem = sd(v) / sqrt(length(v)),
                                      n = length(v)))
stats <- data.frame(timepoint = agg$timepoint, level = agg$sfRNA_level,
                    mean = agg$value[, "mean"], sem = agg$value[, "sem"])
stats$x <- ifelse(stats$timepoint == "24h",
                  match(stats$level, c("L", "M", "H")) - 1,
                  match(stats$level, c("L", "M", "H")) + 2.2)
stats$fill <- COLORS[stats$level]
stats$letter <- LETTERS[paste0(stats$timepoint, ".", stats$level)]
W <- 0.72

# 48h 柱的 45° 斜线纹理(手工计算线段, 裁剪到柱体矩形内)
hatch_segments <- function(x0, x1, y0, y1, spacing = 0.14) {
  segs <- data.frame(x = numeric(0), y = numeric(0),
                     xend = numeric(0), yend = numeric(0))
  for (cc in seq(x0 - y1, x1 - y0, by = spacing)) {
    pts <- rbind(c(x0, x0 - cc), c(x1, x1 - cc),
                 c(cc + y0, y0), c(cc + y1, y1))
    keep <- pts[, 1] >= x0 - 1e-9 & pts[, 1] <= x1 + 1e-9 &
            pts[, 2] >= y0 - 1e-9 & pts[, 2] <= y1 + 1e-9
    pts <- pts[keep, , drop = FALSE]
    if (nrow(pts) >= 2) {
      d <- pts[order(pts[, 1]), , drop = FALSE]
      segs <- rbind(segs, data.frame(x = d[1, 1], y = d[1, 2],
                                     xend = d[nrow(d), 1], yend = d[nrow(d), 2]))
    }
  }
  segs
}
h48 <- do.call(rbind, lapply(which(stats$timepoint == "48h"), function(i) {
  hatch_segments(stats$x[i] - W / 2, stats$x[i] + W / 2, 0, stats$mean[i])
}))

set.seed(42)
pts <- do.call(rbind, lapply(seq_len(nrow(df)), function(i) {
  r <- df[i, ]
  xx <- stats$x[stats$timepoint == r$timepoint & stats$level == r$sfRNA_level]
  data.frame(x = xx + runif(1, -0.14, 0.14), y = r$value)
}))

p <- ggplot() +
  geom_col(data = stats, aes(x = x, y = mean), fill = stats$fill,
           color = "black", linewidth = 0.9, width = W) +
  geom_segment(data = h48, aes(x = x, y = y, xend = xend, yend = yend),
               color = "black", linewidth = 0.55) +
  geom_errorbar(data = stats, aes(x = x, ymin = mean - sem, ymax = mean + sem),
                width = 0.16, linewidth = 0.9, color = "black") +
  geom_point(data = pts, aes(x = x, y = y), shape = 21, size = 3.2,
             fill = "white", color = "black", stroke = 0.9) +
  geom_hline(yintercept = MOCK_VALUE, linetype = "dashed", linewidth = 0.9) +
  geom_text(data = stats, aes(x = x, y = mean + sem + 0.6, label = letter),
            family = "noto", size = 5) +
  annotate("text", x = c(1, 4.2), y = -1.55, label = c("24", "48"),
           family = "noto", size = 5) +
  annotate("text", x = c(1, 4.2), y = -2.9, label = c("hpi", "hpi"),
           family = "noto", size = 5) +
  annotate("text", x = stats$x, y = -0.5, label = stats$level,
           family = "noto", size = 5) +
  annotate("segment", x = -0.7, xend = 5.9, y = 0, yend = 0, linewidth = 0.9) +
  annotate("segment", x = -1.0, xend = -1.0, y = 0, yend = 16.2, linewidth = 0.9) +
  # "a" 与 "Huh7": 单字重 TTF 无粗体, 用 4 次微偏移叠印仿粗
  annotate("text", x = -0.78 + c(-0.008, 0.008, 0, 0), y = 17.3 + c(0, 0, -0.012, 0.012),
           label = "a", family = "noto", size = 9) +
  annotate("text", x = 1.05 + c(-0.008, 0.008, 0, 0), y = 18.0 + c(0, 0, -0.012, 0.012),
           label = "Huh7", family = "noto", size = 8) +
  annotate("text", x = 0.55, y = 13.9, label = "IFN-\u03b2", family = "noto", size = 7) +
  scale_x_continuous(breaks = NULL, limits = c(-1.0, 5.9), expand = c(0, 0)) +
  scale_y_continuous(breaks = c(0, 5, 10, 15), limits = c(-3.6, 18.6),
                     expand = c(0, 0)) +
  labs(y = "Relative gene\nexpression", x = NULL) +
  theme_classic(base_family = "noto", base_size = 13) +
  theme(axis.line = element_line(linewidth = 0.9),
        axis.line.x = element_blank(),
        axis.line.y = element_blank(),
        axis.ticks.x = element_blank(),
        axis.ticks.y = element_line(linewidth = 0.9),
        axis.ticks.length = unit(4, "pt"),
        axis.text = element_text(size = 12, color = "black"),
        axis.text.x = element_blank(),
        axis.title.y = element_text(size = 14, color = "black"))

agg_png(file.path(here, "reproduce_r.png"), width = 6.4, height = 6.6,
        units = "in", res = 300)
print(p)
invisible(dev.off())
cat("saved reproduce_r.png\n")