生化亲和力实验中绘制Sigmoid曲线并获取Kd值的技术问题

一、数据预处理（转换为Tidy格式） #

先将宽格式实验数据转换为适合拟合和可视化的长格式，同时简化单位转换逻辑：

# 加载所需工具包
library(tidyverse)
library(drc)

# 构建示例数据框（可替换为你的read.csv读取逻辑）
mydata <- tibble(
  ligand_conc = c(5.289E-09, 1.058E-08, 2.115E-08, 4.231E-08, 8.462E-08, 1.692E-07, 3.385E-07, 6.769E-07, 1.354E-06, 2.708E-06, 5.415E-06, 1.083E-05, 2.166E-05, 4.332E-05, 8.665E-05, 1.733E-04),
  exp_a = c(792.6444, 792.8537, 788.0273, 793.9693, 792.3848, 792.311, 790.5109, 790.4974, 796.1723, 790.8627, 790.2954, 784.7171, 773.0447, 760.8085, 745.5512, 738.3463),
  exp_b = c(790.2453, 793.8565, 789.5286, 791.8368, 788.5138, 790.0382, 792.85, 789.1439, 790.3487, 792.1872, 786.6738, 780.0627, 775.8658, 762.8376, 747.4288, 737.8717),
  exp_c = c(788.2453, 790.5648, 792.8529, 790.1368, 793.5138, 791.7038, 788.85, 791.1439, 789.4487, 788.8872, 789.5674, 783.3063, 774.8658, 764.5838, 749.4288, 736.5872)
) %>%
  # 转换为长格式，适配ggplot分组逻辑
  pivot_longer(cols = starts_with("exp"), names_to = "experiment", values_to = "Fnorm") %>%
  # 单位转换：M → μM（直接乘以1e6，替代两次*1000）
  mutate(ligand = ligand_conc * 1e6)

二、四参数Sigmoid模型拟合（解决收敛+提取Kd+拟合优度） #

使用drc包的四参数logistic模型（L.4()）拟合数据，该模型的ED50值完全等同于实验Kd值（拐点处的配体浓度）。同时手动计算McFadden's R²评估拟合效果：

# 拟合四参数logistic模型，指定初始值解决收敛问题
sig_model <- drm(Fnorm ~ ligand, data = mydata, fct = L.4(), 
                 start = c(730, 795, 0.5, 0.05)) # 初始值：下限、上限、斜率、ED50

# 提取Kd（ED50）及95%置信区间
kd_summary <- ED(sig_model, c(50), interval = "confidence")
cat("Kd值（95%置信区间）：\n")
print(kd_summary)

# 计算McFadden's R²
null_model <- lm(Fnorm ~ 1, data = mydata)
mcfadden_r2 <- 1 - logLik(sig_model)/logLik(null_model)
cat("\nMcFadden's R²：", round(as.numeric(mcfadden_r2), 4), "\n")

三、生成平滑曲线+95%置信区间的预测数据 #

手动生成覆盖现有浓度范围（含外推区间）的预测数据，规避geom_smooth的置信区间报错问题：

# 生成用于绘制曲线的x值（包含外推的低/高浓度，解决平台缺失问题）
pred_x <- seq(min(mydata$ligand)*0.1, max(mydata$ligand)*10, length.out = 1000)

# 预测y值及95%置信区间
pred_data <- predict(sig_model, newdata = data.frame(ligand = pred_x), 
                     interval = "confidence", level = 0.95) %>%
  as.data.frame() %>%
  mutate(ligand = pred_x) %>%
  rename(Fnorm_fit = Fit, Fnorm_lower = Lower, Fnorm_upper = Upper)

四、ggplot2可视化（平滑曲线+置信区间+数据点） #

替换geom_smooth为手动绘制的曲线和置信区间，解决曲线粗糙、无CI的问题：

ggplot() +
  # 绘制95%置信区间
  geom_ribbon(data = pred_data, aes(x = ligand, ymin = Fnorm_lower, ymax = Fnorm_upper), 
              fill = "grey70", alpha = 0.3) +
  # 绘制平滑拟合曲线
  geom_line(data = pred_data, aes(x = ligand, y = Fnorm_fit), color = "black", linewidth = 1) +
  # 绘制原始实验数据点
  geom_point(data = mydata, aes(x = ligand, y = Fnorm, color = experiment), size = 2) +
  
  # 主题与坐标轴设置
  theme_bw() +
  theme(
    axis.line.x.bottom = element_line(color = 'black'),
    axis.line.y.left = element_line(color = 'black'),
    axis.line.y.right = element_line(color = 'black'),
    panel.grid.minor.x = element_blank(),
    panel.border = element_blank(),
    axis.title.x = element_markdown(),
    axis.title.y = element_markdown(),
    axis.minor.ticks.length = rel(1),
    axis.text = element_text(color = "black", size = 10),
    axis.ticks = element_line(linewidth = 0.6),
    axis.ticks.length = unit(2.75, "pt")
  ) +
  
  # 分组颜色设置
  scale_color_manual(
    name = "Replicate",
    labels = c("1", "2", "3"),
    values = c("sienna1", "dodgerblue", "grey43")
  ) +
  
  # 坐标轴范围与转换
  coord_cartesian(ylim = c(720, 800), expand = TRUE) +
  scale_x_continuous(
    trans = "log10",
    expand = c(0, 0),
    label = label_number(),
    breaks = c(0.01, 0.1, 1, 10, 100, 1000),
    guide = guide_axis_logticks(long = 2.3, mid = 1.65, short = 0.75),
    limits = c(0.001, 1000)
  ) +
  
  # 标签设置
  labs(
    y = "Fnorm [%<sub>280</sub>]",
    x = "Ligand conc. [µM]"
  )

关键问题说明 #

1. 曲线粗糙：通过生成1000个预测点替代原始数据点，确保曲线平滑。
2. 置信区间报错：geom_smooth对drm的CI支持不完善，改用predict手动计算后绘制。
3. ED50与Kd的关系：四参数模型的ED50是响应值达到中间水平对应的浓度，完全匹配实验Kd定义。
4. 平台缺失外推：通过pred_x扩展浓度范围，模型会基于拟合的上下限参数生成完整Sigmoid曲线。
5. 收敛失败：指定start参数给出合理初始值（上下限、斜率、ED50），帮助模型找到最优解。

内容的提问来源于stack exchange，提问作者Cobalamin