通过原点的多项式回归线与校准曲线中的方程

Question

我想要一条通过零绘制的二阶（是吗）回归线，关键是我需要这种关系的方程。

这是我的数据：

ecoli_ug_ml  A420 rpt
1            0 0.000   1
2           10 0.129   1
3           20 0.257   1
4           30 0.379   1
5           40 0.479   1
6           50 0.579   1
7           60 0.673   1
8           70 0.758   1
9           80 0.838   1
10          90 0.912   1
11         100 0.976   1
12           0 0.000   2
13          10 0.126   2
14          20 0.257   2
15          30 0.382   2
16          40 0.490   2
17          50 0.592   2
18          60 0.684   2
19          70 0.772   2
20          80 0.847   2
21          90 0.917   2
22         100 0.977   2
23           0 0.000   3
24          10 0.125   3
25          20 0.258   3
26          30 0.376   3
27          40 0.488   3
28          50 0.582   3
29          60 0.681   3
30          70 0.768   3
31          80 0.846   3
32          90 0.915   3
33         100 0.977   3

我的情节是这样的：（sci2 只是一些轴和文本格式，必要时可以包含）

ggplot(calib, aes(ecoli_ug_ml, A420)) +
geom_point(shape=calib$rpt) +
stat_smooth(method="lm", formula=y~poly(x - 1,2)) +
scale_x_continuous(expression(paste(italic("E. coli"),~"concentration, " ,mu,g~mL^-1,))) +
scale_y_continuous(expression(paste(Absorbance["420nm"], ~ ", a.u."))) +
sci2

当我看到这个时，这条线与点的拟合非常好。

当我检查 coef 时，我 认为 y 截距不为零（这对我的目的来说是不可接受的）但老实说我并不真正理解这些行：

coef(lm(A420 ~ poly(ecoli_ug_ml, 2, raw=TRUE), data = calib))                 
                  (Intercept) poly(ecoli_ug_ml, 2, raw = TRUE)1 
                -1.979021e-03                      1.374789e-02 
poly(ecoli_ug_ml, 2, raw = TRUE)2 
                -3.964258e-05 

因此，我认为情节实际上也不完全正确。

所以，我需要生成一条强制通过零的回归线并得到它的方程，并且，当它给我所说的方程时，理解它在说什么。如果我可以直接用方程式注释绘图区域，我会非常激动。

我现在花了大约 8 个小时来解决这个问题，我检查了 excel 并在 8 秒内得到了一个公式，但我真的很想使用 R 来解决这个问题。谢谢！

澄清一下：这个图的主要目的不是展示这些数据的分布，而是提供一个视觉确认，即我从这些点生成的方程非常适合读数

Answer 1

summary(lm(A420~poly(ecoli_ug_ml,2,raw=T),data=calib))

# Call:
# lm(formula = A420 ~ poly(ecoli_ug_ml, 2, raw = T), data = calib)
# ...
# Coefficients:
#                                  Estimate Std. Error t value Pr(>|t|)    
# (Intercept)                    -1.979e-03  1.926e-03  -1.028    0.312    
# poly(ecoli_ug_ml, 2, raw = T)1  1.375e-02  8.961e-05 153.419   <2e-16 ***
# poly(ecoli_ug_ml, 2, raw = T)2 -3.964e-05  8.631e-07 -45.932   <2e-16 ***
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

# Residual standard error: 0.004379 on 30 degrees of freedom
# Multiple R-squared:  0.9998,  Adjusted R-squared:  0.9998 
# F-statistic: 8.343e+04 on 2 and 30 DF,  p-value: < 2.2e-16

因此截距不完全为 0，但与 Std 相比它很小。错误。也就是说，截距与0没有显着差异。

您可以通过这种方式强制拟合而不使用截距（注意公式中的 -1）：

summary(lm(A420~poly(ecoli_ug_ml,2,raw=T)-1,data=calib))

# Call:
# lm(formula = A420 ~ poly(ecoli_ug_ml, 2, raw = T) - 1, data = calib)
# ...
# Coefficients:
#                                  Estimate Std. Error t value Pr(>|t|)    
# poly(ecoli_ug_ml, 2, raw = T)1  1.367e-02  5.188e-05  263.54   <2e-16 ***
# poly(ecoli_ug_ml, 2, raw = T)2 -3.905e-05  6.396e-07  -61.05   <2e-16 ***
# ---
# Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

# Residual standard error: 0.004383 on 31 degrees of freedom
# Multiple R-squared:      1,   Adjusted R-squared:      1 
# F-statistic: 3.4e+05 on 2 and 31 DF,  p-value: < 2.2e-16

请注意，系数没有明显变化。

编辑（回应 OP 的评论）

stat_smooth(...) 中指定的公式只是直接传递给lm(...) 函数，因此您可以在stat_smooth(...) 中指定任何适用于lm(...) 的公式。上述结果的要点是，即使不强制截距为 0，与 y (0-1) 中的范围相比，它也非常小 (-2e-3)，因此绘制有和没有的曲线将给出几乎无法区分的结果.您可以通过运行以下代码自己查看：

ggplot(calib, aes(ecoli_ug_ml, A420)) +
  geom_point(shape=calib$rpt) +
  stat_smooth(method="lm", formula=y~poly(x,2,raw=T),colour="red") +
  stat_smooth(method="lm", formula=y~-1+poly(x,2,raw=T),colour="blue") +
  scale_x_continuous(expression(paste(italic("E. coli"),~"concentration, " ,mu,g~mL^-1,))) +
  scale_y_continuous(expression(paste(Absorbance["420nm"], ~ ", a.u.")))

蓝色和红色曲线几乎相互重叠，但并不完全重叠（您可能需要打开绘图窗口才能看到它）。不，您确实不必须在“ggplot 之外”执行此操作。

您报告的问题与使用默认 raw=F 有关。这导致poly(...) 使用正交多项式，根据定义具有常数项。所以使用y~-1+poly(x,2) 确实没有意义，而使用y~-1+poly(x,2,raw=T) 确实有意义。

最后，如果所有使用poly(...) 或不使用raw=T 的业务都会引起混淆，那么使用formula = y~ -1 + x + I(x^2) 可以达到完全相同的结果。这适合二阶多项式 (a*x +b*x^2) 并抑制常数项。

Answer 2

我认为您误解了 Intercept 以及 stat_smooth 的工作原理。统计学家完成的多项式拟合通常不使用 raw=TRUE 参数。默认值为 FALSE，多项式被构造为正交的，以便在查看标准误差时对拟合改进进行适当的统计评估。如果您尝试通过在公式中使用-1 或0+ 来消除拦截，看看会发生什么是很有启发性的。尝试使用您的数据和代码来摆脱拦截：

   ....+
   stat_smooth(method="lm", formula=y~0+poly(x - 1,2)) + ...

您将看到拟合线在 -0.5 处与 y 轴相交并发生变化。现在看看截距的非原始值。

  coef(lm(A420~poly(ecoli_ug_ml,2),data=ecoli))
      (Intercept) poly(ecoli_ug_ml, 2)1 poly(ecoli_ug_ml, 2)2 
        0.5466667             1.7772858            -0.2011251 

因此，截距将整条曲线向上移动，以使多项式拟合具有最佳拟合曲率。如果你想用 ggplot2 画一条符合不同规范的线，你应该在 ggplot2 之外计算它，然后在没有误差带的情况下绘制它，因为它确实不具有适当的统计属性。

尽管如此，这是应用在这种情况下是微不足道的调整的方法，我提供它只是为了说明如何添加一组外部派生值。我认为像这样的 _ad_hoc_ 调整在实践中是危险的：

 mod <-   lm(A420~poly(ecoli_ug_ml,2), data=ecoli)
 ecoli$shifted_pred <- predict(mod) - predict( mod, newdata=list(ecoli_ug_ml=0))
 ggplot(ecoli, aes(ecoli_ug_ml, A420)) +
    geom_point(shape=ecoli$rpt)  +
    scale_x_continuous(expression(paste(italic("E. coli"),~"concentration, " ,mu,g~mL^-1,))) +
    scale_y_continuous(expression(paste(Absorbance["420nm"], ~ ", a.u.")))+
    geom_line(data=ecoli, aes(x= ecoli_ug_ml, y=shifted_pred ) )