ggplot 指数平滑，在 exp 内调整参数

Question

ggplot 提供了各种确定趋势线形式的“平滑方法”或“公式”。但是，我不清楚如何指定公式的参数以及如何让指数公式适合我的数据。换句话说，如何告诉 ggplot 它应该适合 exp 中的参数。

df <- data.frame(x = c(65,53,41,32,28,26,23,19))
df$y <- c(4,3,2,8,12,8,20,15)

   x  y
1 65  4
2 53  3
3 41  2
4 32  8
5 28 12
6 26  8
7 23 20
8 19 15
p <- ggplot(data = df, aes(x = x, y = y)) +
  geom_smooth(method = "glm", se=FALSE, color="black", formula = y ~ exp(x)) +
  geom_point()

p

有问题的配合：

但是如果指数内的参数是合适的，那么趋势线的形式就变得合理了：

p <- ggplot(data = df, aes(x = x, y = y)) +
  geom_smooth(method = "glm", se=FALSE, color="black", formula = y ~ exp(-0.09 * x)) +
  geom_point()

p

Answer 1

这是一种使用方法nls 而不是glm 的方法。

您可以使用method.args = 中提供的列表将其他参数传递给nls。在这里，我们定义了要拟合的 a 和 r 系数的起始值。

library(ggplot2)
ggplot(data = df, aes(x = x, y = y)) +
  geom_smooth(method = "nls", se = FALSE,
              formula = y ~ a * exp(r * x),
              method.args = list(start = c(a = 10, r = -0.01)),
              color = "black") +
  geom_point()

正如 cmets 中所讨论的，在图表上获取系数的最佳方法是在 ggplot 调用之外拟合模型。

model.coeff <- coef(nls( y ~ a * exp(r * x), data = df, start = c(a = 50, r = -0.04)))

ggplot(data = df, aes(x = x, y = y)) +
  geom_smooth(method = "nls", se = FALSE,
              formula = y ~ a * exp(r * x),
              method.args = list(start = c(a = 50, r = -0.04)),
              color = "black") +
  geom_point() + 
  geom_text(x = 40, y = 15,
            label = as.expression(substitute(italic(y) == a %.% italic(e)^(r %.% x),
                                             list(a = format(unname(model.coeff["a"]),digits = 3),
                                                  r = format(unname(model.coeff["r"]),digits = 3)))),
            parse = TRUE)

Answer 2

首先，要将附加参数传递给传递给geom_smooth 的method 参数的函数，您可以将命名参数列表传递给method.args。

其次，您看到的问题是 glm 将系数放在整个术语的前面：y ~ coef * exp(x) 而不是里面：y ~ exp(coef * x) 就像你想要的那样。您可以使用优化来解决 glm 之外的后者，但您可以通过转换将其放入 GLM 范式：日志链接。这是有效的，因为它就像采用您想要拟合的方程 y = exp(coef * x) 并获取双方的日志，所以您现在正在拟合 log(y) = coef * x，这相当于您想要拟合并使用 GLM 范式. （这忽略了截距。它也以转换后的链接单元结束，但如果你愿意，它很容易转换回来。）

您可以在 ggplot 之外运行它以查看模型的外观：

df <- data.frame(
    x = c(65,53,41,32,28,26,23,19), 
    y <- c(4,3,2,8,12,8,20,15)
)

bad_model <- glm(y ~ exp(x), family = gaussian(link = 'identity'), data = df)
good_model <- glm(y ~ x, family = gaussian(link = 'log'), data = df)

# this is bad
summary(bad_model)
#> 
#> Call:
#> glm(formula = y ~ exp(x), family = gaussian(link = "identity"), 
#>     data = df)
#> 
#> Deviance Residuals: 
#>     Min       1Q   Median       3Q      Max  
#> -7.7143  -2.9643  -0.8571   3.0357  10.2857  
#> 
#> Coefficients:
#>               Estimate Std. Error t value Pr(>|t|)   
#> (Intercept)  9.714e+00  2.437e+00   3.986  0.00723 **
#> exp(x)      -3.372e-28  4.067e-28  -0.829  0.43881   
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> (Dispersion parameter for gaussian family taken to be 41.57135)
#> 
#>     Null deviance: 278.00  on 7  degrees of freedom
#> Residual deviance: 249.43  on 6  degrees of freedom
#> AIC: 56.221
#> 
#> Number of Fisher Scoring iterations: 2

# this is better
summary(good_model)
#> 
#> Call:
#> glm(formula = y ~ x, family = gaussian(link = "log"), data = df)
#> 
#> Deviance Residuals: 
#>    Min      1Q  Median      3Q     Max  
#> -3.745  -2.600   0.046   1.812   6.080  
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)  3.93579    0.51361   7.663 0.000258 ***
#> x           -0.05663    0.02054  -2.757 0.032997 *  
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> (Dispersion parameter for gaussian family taken to be 12.6906)
#> 
#>     Null deviance: 278.000  on 7  degrees of freedom
#> Residual deviance:  76.143  on 6  degrees of freedom
#> AIC: 46.728
#> 
#> Number of Fisher Scoring iterations: 6

从这里，您可以重现 geom_smooth 将要做的事情：在整个域中创建一系列 x 值，并将预测用作该行的 y 值：

# new data is a sequence across the domain of the model
new_df <- data.frame(x = seq(min(df$x), max(df$x), length = 501))

# `type = 'response'` because we want values for y back in y units
new_df$bad_pred <- predict(bad_model, newdata = new_df, type = 'response')
new_df$good_pred <- predict(good_model, newdata = new_df, type = 'response')

library(tidyr)
library(ggplot2)

new_df %>% 
    # reshape to long form for ggplot
    gather(model, y, contains('pred')) %>%
    ggplot(aes(x, y)) + 
    geom_line(aes(color = model)) + 
    # plot original points on top
    geom_point(data = df)

当然，让 ggplot 为您处理这一切要容易得多：

ggplot(df, aes(x, y)) + 
    geom_smooth(
        method = 'glm', 
        formula = y ~ x, 
        method.args = list(family = gaussian(link = 'log'))
    ) + 
    geom_point()