# Bland–Altman 一致性图：两种测量方法的差值与一致性界限
# 项目原创脚本，MIT 许可；模拟数据 CC0-1.0；图形 CC BY 4.0。
# 在全新 R 会话运行。每次运行会更新下方 output_dir 中的本示例文件。
# 图形标签使用英文以避免不同电脑缺少中文字体；可在参数和 labs() 中修改。

# 1. 依赖：首次使用时手动运行以下安装命令，脚本不会自动安装。
# install.packages(c("ggplot2", "jsonlite"))
required_packages <- c("ggplot2", "jsonlite")
missing_packages <- required_packages[!vapply(required_packages, requireNamespace, logical(1), quietly = TRUE)]
if (length(missing_packages)) stop("Please install: ", paste(missing_packages, collapse = ", "))
library(ggplot2)

# 2. 可修改参数：颜色、字体、导出尺寸和路径。输出路径相对于当前工作目录。
set.seed(20261014)
colors <- c("#007F96", "#B66A46")
font_family <- "sans"
figure_width <- 8
figure_height <- 5.5
output_dir <- file.path(getwd(), "rfig-bland-altman-output")
dir.create(output_dir, recursive = TRUE, showWarnings = FALSE)
theme_set(theme_minimal(base_size = 12, base_family = font_family) +
 theme(plot.background = element_rect(fill = "white", color = NA), panel.grid.minor = element_blank(),
       plot.title = element_text(face = "bold", size = 15), legend.position = "bottom",
       plot.margin = margin(14, 18, 14, 14)))

# 3. 示例数据：每行一个对象，两种方法各测量一次。 保持两方法按受试者正确配对，检查单位相同，并预先确定临床可接受差值。
# 使用自己的 CSV 时，将以下数据生成段替换为：
# data <- read.csv("your-data.csv", stringsAsFactors = FALSE)
# 再按本例要求设置分类顺序、检查字段、单位及缺失值。
truth <- rnorm(70, 50, 8)
data <- data.frame(id = 1:70, method_a = truth + rnorm(70, 0, 2), method_b = truth + 1.5 + rnorm(70, 0, 2))
stopifnot(!anyNA(data))

# 4. 统计计算与绘图。差值定义为 B−A；横轴为两方法均值。偏倚为平均差，一致性界限为平均差 ± 1.96×差值 SD，近似覆盖 95% 差值；不是偏倚或界限的 95% CI。需检查差值近似正态、异方差及比例偏倚；可接受范围由专业背景预先确定。
data$average <- (data$method_a + data$method_b) / 2
data$difference <- data$method_b - data$method_a
bias <- mean(data$difference); lower <- bias - 1.96 * sd(data$difference); upper <- bias + 1.96 * sd(data$difference)
lines <- data.frame(y = c(lower, bias, upper), label = c(sprintf("Lower LoA %.2f", lower), sprintf("Bias %.2f", bias), sprintf("Upper LoA %.2f", upper)), type = c("LoA", "Bias", "LoA"))
p <- ggplot(data, aes(average, difference)) + geom_point(color = colors[1], alpha = .7, size = 2.2) +
  geom_hline(data = lines, aes(yintercept = y, linetype = type), color = colors[2]) +
  geom_text(data = lines, aes(x = Inf, y = y, label = label), inherit.aes = FALSE, hjust = 1.05, vjust = -0.5, size = 3.3) +
  scale_linetype_manual(values = c(Bias = "solid", LoA = "dashed")) + guides(linetype = "none") +
  scale_y_continuous(expand = expansion(mult = c(.12, .18))) + labs(x = "Mean of methods A and B (demo units)", y = "Difference B - A (demo units)", title = "Agreement between two measurement methods")
metrics <- list(bias = bias, lower = lower, upper = upper)
print(p)

# 5. 导出：CSV 与绘图所用 data 完全一致，preview 与科研图片均来自同一个 p。
write.csv(data, file.path(output_dir, "data.csv"), row.names = FALSE, fileEncoding = "UTF-8")
ggsave(file.path(output_dir, "preview.png"), p, width = figure_width, height = figure_height, dpi = 160, bg = "white")
ggsave(file.path(output_dir, "figure.png"), p, width = figure_width, height = figure_height, dpi = 300, bg = "white")
ggsave(file.path(output_dir, "figure.pdf"), p, width = figure_width, height = figure_height, device = grDevices::pdf, bg = "white")
ggsave(file.path(output_dir, "figure.tiff"), p, width = figure_width, height = figure_height, dpi = 300, compression = "lzw", bg = "white")
jsonlite::write_json(metrics, file.path(output_dir, "metrics.json"), auto_unbox = TRUE, digits = 16, pretty = TRUE)